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Abstract 


Reachability  and  Real-Time  Actuation  Strategies  for  the  Active  SLIP  Model 

by 

Giulia  Piovan 


Running  and  hopping  follow  similar  patterns  for  different  animals,  independent  of  the 
number  of  legs  employed.  An  aerial  phase  alternates  with  a  ground  contact  phase,  dur¬ 
ing  which  the  center  of  mass  moves  as  if  a  spring  were  compressed  and  then  extended 
to  recover  stored  elastic  energy.  Hence,  consisting  of  a  point  mass  mounted  on  a  mass¬ 
less  spring  leg,  the  Spring  Loaded  Inverted  Pendulum  (SLIP)  is  a  prevalent  model  for 
analyzing  running  and  hopping.  In  this  work  we  consider  an  actuated  version  of  the 
SLIP  model,  with  a  series  elastic  actuator  added  to  the  leg,  serving  the  purposes  of 
adding/removing  energy  to/from  the  system  and  of  modifying  dynamics  during  stance, 
toward  achieving  non-steady  locomotion  on  varying  terrain.  While  the  SLIP  model  has 
been  a  topic  of  research  in  legged  locomotion  for  several  decades,  studies  on  the  effect 
of  actuation  on  the  system’s  behavior  are  still  not  complete.  The  goal  of  this  thesis 
is  to  explore  how  a  series  elastic  actuator  applied  to  the  SLIP  model’s  leg  can  change 
the  system’s  dynamics.  This,  in  turn,  enables  a  variety  of  long-term  planning  strategies 
for  using  limited  footholds  and  design  non-steady  gaits  while  simultaneously  recovering 
from  unexpected  perturbations,  both  sensorial  and  due  to  a  limited  knowledge  of  the 
terrain  profile.  We  principally  investigate  how,  through  actuation,  we  can  solve  partially 
or  completely  the  system’s  equations  of  motion,  to  enforce  a  desired  trajectory  and  reach 
a  desired  state.  We  also  determine  the  reachable  state  space  of  the  model  using  several 
different  actuation  strategies,  investigating  the  variation  of  the  reachable  set  with  respect 

to  particular  actuator  motions  and  providing  relationships  between  local  actuator  dis- 
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placements  throughout  stance  and  location  of  the  reached  apex  state.  We  then  propose 
a  control  strategy  based  on  graphical  and  numerical  studies  of  the  reachability  space  to 
drive  the  system  to  a  desired  state,  with  the  ability  to  reduce  the  effects  of  sensing  errors 
and  disturbances  happening  at  landing  as  well  as  during  ground  contact. 
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Chapter  1 


Introduction  and  Motivation 


Thanks  to  consolidated  knowledge  in  modelling,  manufacturing  and  control,  wheeled 
vehicles  are  still  the  most  popular  choice  in  terms  of  locomotion.  Their  efficiency  on 
continuous  surfaces,  however,  is  paired  with  a  lack  of  versatility  and  adaptability  that 
poses  a  limitation  in  many  real-world  applications,  particularly  those  instances  that  see 
the  presence  of  obstacles  in  the  terrains.  These  same  terrains  can  easily  be  reached 
by  most  animals  and  humans.  Hence,  legged  locomotion  has  been  a  topic  of  study 
in  robotics  and  has  seen  remarkable  advances  in  the  last  forty  years.  The  advantages  of 
legged  systems  over  wheeled  or  track-based  systems  are  multifold.  Firstly,  the  continuous 
path  of  support  needed  by  wheeled  vehicles  is  substituted  by  a  set  of  isolated  footholds  in 
the  case  of  legged  systems  and,  as  a  consequence,  not  all  the  terrain  needs  to  be  specified 
and  to  be  accessible.  This  implies  an  ability  to  adapt  to  surfaces  that  are  uneven  or  with 
big  variations.  Additionally,  the  movement  is  determined  by  the  robot’s  dynamics  more 
than  by  the  terrain  surface  due  to  an  active  suspension  of  the  body.  Finally,  by  studying 
legged  locomotion,  we  hope  to  gain  insights  into  human  mobility  as  well,  which  can  be 
applied  to  different  areas  of  medicine,  for  example  prosthetics. 

Based  on  our  daily  experience,  walking  on  legs  seems  to  be  a  trivial  task,  and  yet  the 
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capabilities  of  legged  animals  are  far  superior  to  those  of  any  existing  robot  in  that  animals 
can  negotiate  a  variety  of  different  terrains  in  an  extremely  efficient  way.  However, 
since  the  mechanisms  of  locomotion  are  yet  not  fully  understood,  and  are  coupled  with 
mechanical  and  sensorial  limitations  ([1]),  legged  locomotion  is  still  very  much  an  open 
problem.  First  of  all,  it  is  governed  by  a  much  more  complicated  dynamics  than  wheeled 
systems,  which  presents  a  challenge  from  the  controller’s  point  of  view.  Furthermore, 
legged  systems  have  higher  peak  energy  and  torque  requirements  than  wheeled  vehicles, 
and  they  experience  energy  loss  at  each  impact  and  due  to  friction  at  the  joints.  On 
hardware  platforms,  autonomous  legged  robots  need  to  carry  all  their  weight  (motors, 
actuators,  batteries)  on  their  legs,  which  have  to  be  robust. 

1.1  Literature  Review 

The  two  basic  gaits  available  to  legged  systems  are  walking  and  running.  Loosely 
speaking,  the  distinction  between  the  two  is  dependent  on  the  presence  or  absence  of  an 
aerial  phase:  while  running,  there  is  at  most  one  leg  in  contact  with  the  ground  at  any 
given  time,  followed  by  a  state  with  no  ground  contact,  Fig.  1.1. 


Figure  1.1:  Man  running  at  full  speed.  From  Muybridge’s  Human  Figure  in  Motion 
(1872-1885),  by  Eadweard  Muybridge.  Courtesy  of  Dover  Publications,  Inc. 

Legs  acting  in  unison  in  some  respects  can  be  considered  as  equivalent  to  a  virtual 
leg  located  at  the  center  of  the  pair.  For  this  reason,  numerous  researchers  interested  in 
running  have  focused  their  attention  to  single-legged  systems. 
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The  conventional  template  used  for  analyzing  running  and  hopping  is  the  so-called 
Spring  Loaded  Inverted  Pendulum  (SLIP)  model.  Consisting  of  a  point  mass  mounted  on 
a  massless  spring  leg,  this  2  degrees  of  freedom  (DOF)  system  is  energetically  conserva¬ 
tive.  The  2-DOF  SLIP  model  can  be  extended  to  its  3-DOF  counterpart  by  substituting 
the  point  mass  with  a  body  that  is  free  to  rotate  via  an  applied  torque  at  the  hip.  The 
spring-mass  system  mimics  animals  that  hop  at  a  particular  preferred  frequency  below 
which  the  motion  requires  more  energy.  Furthermore,  it  can  be  seen  that  hopping  fol¬ 
lows  similar  patterns  for  different  animals,  e.g.,  the  ratio  of  contact  time  versus  aerial 
time,  which  are  independent  on  the  number  of  legs.  During  contact  with  the  ground, 
the  center  of  mass  moves  as  if  a  spring  were  compressed  and  then  extended  recovering 
its  stored  elastic  energy.  Hence,  in  combination  with  biological  data  fitting,  this  simple 
model  can  predict  the  center  of  mass  dynamics  of  running  insects  and  animals,  as  shown 
in  [2],  in  [3],  as  well  in  [4],  [5],  and  [6].  These  works  stress  the  capability  of  the  SLIP 
model  of  describing  the  interdependency  between  physical  and  morphological  parameters 
characterizing  running  and  hopping,  such  as  frequency,  stiffness  and  stride  time.  In  gen¬ 
eral,  animals’  anatomy  is  far  more  complex  than  what  the  SLIP  model  captures,  due  to 
the  presence  of  joints,  ankles,  knees,  hips,  leg  mass  and  inertia,  and  muscles  that  act  as 
nonlinear  springs.  However,  as  Full  and  Koditschek  explain  in  [1],  the  SLIP  model  is  a 
valid  template  which  removes  the  non  essential  complexity  of  the  animal’s  anatomy  but 
still  preserves  the  behavior  that  serves  as  a  guide  for  the  development  of  control  strategies 
for  locomotion.  These  strategies  can  then  be  applied  to  an  anchor,  a  more  complex  and 
realistic  model  to  mimic  the  animal’s  morphology  and  physiology. 

Studies  on  guinea  fowls  ([7]),  and  on  insects  ([8],  and  [9]),  show  that  animals  and 
insects  use  their  legs’  natural  actuation  to  stabilize  their  body.  This  suggests  that  the 
passive  SLIP  is  not  fully  able  to  model  running  and  hopping  gaits.  Since  it  is  clearly 
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desirable  to  enable  legged  locomotion  in  situations  with  variability  in  both  terrain  height 
and  forward  velocity,  requiring  accompanying  modification  of  net  energy,  this  justifies 
the  introduction  of  an  actuated  version  of  the  SLIP.  There  are  several  methods  to  add 
actuation  to  the  SLIP  model.  The  most  common  is  the  use  of  a  telescopic  leg  with  elastic 
actuator  acting  during  the  stance  phase  to  modify  the  net  energy  of  the  system  ([10], 
[11],  [12], [13],  [14],  [15],  [16]).  Authors  of  [17]  added  a  clock-based  torque  at  the  hip  of 
the  system,  while  [18]  propose  instead  to  modify  the  leg  length  during  the  flight  phase. 
In  [19],  authors  consider  the  effect  of  modifying  both  leg  length  and  spring  stiffness  during 
stance,  limiting  their  analysis  to  the  vertical  hopping  case.  These  added  actuations  have 
been  used  to  reach  different  energy  levels,  meaning  changing  the  velocity  or  height  of  the 
body,  toward  controlling  step  size. 

1.1.1  Control  Strategies 

The  passive  2-DOF  SLIP  model  has  been  employed  as  a  starting  point  in  studies 
focusing  on  control  strategies  to  stabilize  the  system  in  the  presence  of  disturbances  due 
to  uneven  terrain  or  unexpected  forces.  Because  of  the  completely  passive  nature  of  the 
system,  its  main  tool  in  terms  of  balance  is  the  placement  of  the  swing  leg.  In  general, 
to  each  forward  speed  corresponds  a  ground  contact  leg  angle  for  steady-state  motion. 
Variations  of  this  angle  have  the  effect  of  decelerating  or  accelerating  the  system.  A 
proper  adjustment  of  spring  stiffness  and  fixed  leg  angle  serves  as  a  control  strategy  for 
stabilization  [20],  and  in  [21]  it  is  shown  how  choosing  a  fixed  leg  angle  at  touch-down 
results  iu  asymptotically  stable  periodic  gaits. 

By  observing  that  running  animals  retract  their  legs  right  before  contact  to  the 
ground,  Seyfarth  et  al.  ([22])  propose  a  swing- leg  retraction  strategy  for  steady-state 
running.  Instead  of  being  kept  fixed  during  the  flight  phase,  the  leg  angle  moves  as  a 
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linear  function  of  time.  Due  to  the  constant  rotational  speed  at  flight,  the  leg  angle  at 
ground  contact  varies  as  a  function  of  the  apex  state  (i.e.,  the  highest  point  reached  dur¬ 
ing  flight),  and  the  relative  speed  of  the  foot  with  respect  to  the  ground  is  reduced.  The 
system  is  stabilized  at  steady-state  faster  than  when  using  a  fixed  angle  of  attack,  and 
it  is  stable  with  respect  to  variations  on  the  spring  stiffness,  the  forward  speed  and  the 
leg  angle.  The  optimal  retraction  rates  for  maximum  disturbance  rejection,  for  minimal 
energy  loss  at  impact,  and  for  minimal  slipping,  are  studied  by  Karssen  et  al.  in  [23]. 
Return  maps  for  both  fixed  contact  angle  and  swing-leg  retraction  are  computed  in  [24] . 
This  stability  analysis  determines  necessary  parametric  condition  for  asymptotic  stability, 
and  it  is  used  to  quantify  sensory  costs. 

Another  parameter  that  could  be  varied  is  the  leg  stiffness.  Stance  time,  and  conse¬ 
quently  leg  compression,  is  a  function  of  leg  stiffness,  and  therefore  a  careful  choice  can 
maximize  the  passive  stability  of  the  system.  For  example,  in  [25],  Ernst  and  colleagues 
propose  a  feedforward  controller  that  continuously  adapts  spring  stiffness  and  leg  angle 
during  flight  (following  a  precomputed  curve  stored  in  a  look-up  table)  regardless  of  when 
the  foot  touches  the  ground,  thus  allowing  the  system  to  hop  on  uneven  terrain.  Riese 
and  Seyfarth  [19]  investigate  stability  of  vertical  hopping  while  varying  spring  stiffness 
and  leg  length. 

For  the  active  SLIP  model,  the  most  famous  control  strategy  is  proposed  by  Raib- 
ert  [10]:  a  constant  thrust  is  applied  during  ground  contact  at  the  point  of  maximum 
spring  compression  in  order  to  change  the  net  energy  of  the  system.  Theoretical  work 
in  [12],  integrated  by  experimental  work  in  [26],  considers  a  3-DOF  SLIP  model  with 
a  series  elastic  actuator  at  the  leg,  and  a  rotational  actuator  and  a  spring  at  the  hip. 
The  control  action  proposed  exploits  the  system’s  passive  dynamics  motion  in  the  hip  to 
minimize  energy  consumption  and  reach  a  steady  state  gait  at  a  target  forward  speed. 

Other  studies,  either  with  the  passive  or  active  SLIP  model,  focus  on  negotiating 
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uneven  terrains:  e.g.,  [27]  propose  three  strategies  for  step-length  adjustment  that  involve 
changing  either  one  of  forward  speed,  stance  time  or  flight  time,  while  keeping  the  other 
two  constant,  [25]  propose  a  control  action  to  keep  the  running  speed  constant.  Koepl  and 
Hurst  [28]  apply  the  leg  controller  proposed  in  [25]  to  the  leg  actuated  3-DOF  SLIP  model, 
and  incorporate  it  with  a  controller  of  the  impulses  applied  by  leg  and  hip  actuators  during 
ground  contact  to  match  the  impulse  profile  of  the  passive  SLIP  model,  in  order  to  reject 
unexpected  ground  disturbances  (of  damping,  stiffness  and  surface  height).  Schmitt  and 
Clark  [13]  propose  a  leg-length  actuation  strategy  inspired  by  muscle  activations  seen  in 
cockroaches  running  over  uneven  terrain.  The  actuator  is  moved  following  a  sinusoidal 
function  of  time,  so  that  energy  is  removed  from  the  system  during  leg  compression 
and  added  during  leg  extension.  This  latter  control  strategy  is  extended  by  Andrews  et 
al.  [14],  who  compare  its  stability  to  Raibert’s  fixed  thrust  controller  and  to  the  swing- leg 
retraction  technique,  both  in  simulation  and  experimentally. 

The  authors  of  [29]  propose  an  algorithm  for  trajectory  planning,  robust  to  model 
uncertainty  and  measurement  noise. 

Work  in  [15]  combines  an  actuator  displacement  strategy  with  a  model  predictive 
control  approach  to  plan  a  trajectory  based  on  a  fixed  set  of  desired  footholds. 

1.1.2  Approximation  of  the  system’s  dynamics 

While  the  SLIP  model  is  a  very  simple  system,  there  is  not  an  analytic  solution 
for  its  dynamics  during  ground  contact  [30].  Trajectories  have  to  be  pre-computed  and 
stored  in  look-up  tables,  which  can  be  timely  and  computationally  expensive  to  generate, 
store  and  access.  Furthermore,  motions  not  pre-computed  need  to  be  interpolated  from 
existing  trajectories.  Additionally,  the  lack  of  closed-form  solution  represents  a  limitation 
on  our  ability  to  understand  the  SLIP  dynamics  and  to  design  and  implement  real-time 
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controllers,  as  well  as  when  generating  return  maps  to  investigate  the  system’s  stability. 

To  overcome  these  issues,  several  approximations  to  the  closed-form  solution  have 
been  given  throughout  the  years.  Schwind  and  Koditschek  [31]  apply  the  mean  value 
theorem  to  the  integral  of  the  SLIP  dynamics,  and  find  a  closed-form  solution  for  the 
case  where  there  is  conservation  of  momentum.  Then,  the  effects  of  gravity  are  iteratively 
added  to  equations  so  that  under  certain  assumptions  the  solution  will  converge  to  the 
exact  SLIP  dynamics.  Geyer  et  al.  [32]  have  proposed  a  solution  based  upon  conservation 
of  momentum  and  the  hypothesis  of  small  angle  span.  Conservation  of  momentum  is  a 
consequence  of  ignoring  the  effect  of  gravity  along  the  leg;  hence,  the  accuracy  of  this 
solution  degrades  the  more  the  trajectory  deviates  from  the  symmetric  case. 

Locomotion  on  uneven  terrain  may  imply  non-symmetric  strides,  where  the  effect 
of  gravity  cannot  be  neglected.  Following  this  observation,  works  by  Arslan  et  al.  [33] 
and  Yu  et  al.  [34]  propose  approximations  which  include  gravity  corrections,  improving 
the  model  performance.  More  recently,  [35]  propose  an  approximation  of  the  stance 
phase  that  yields  linear  equations,  thus  allowing  to  write  the  entire  SLIP  dynamics  as  a 
piecewise  linear  function. 

1.1.3  SLIP-based  robots 

In  the  80s,  Mark  Raibert  [10]  presented  the  first  and  best-known  planar  monopod 
hopping  robot,  the  so-called  Raibert  hopper,  Fig.  1.2a.  The  Raibert  hopper  consists  of  a 
body  and  a  compliant  elastic  leg,  both  with  mass  and  inertia,  connected  by  a  joint  at  the 
hip.  An  actuator  at  the  hip  can  provide  torque  between  the  body  and  the  leg,  while  a 
series  elastic  actuator  at  the  leg  compresses  or  extends  the  spring.  The  three-part  control 
strategy  considers  hopping  height,  speed  and  body  attitude  as  three  separate  problems: 
forward  velocity  is  controlled  by  the  positioning  of  the  leg;  hopping  height  is  regulated 
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by  the  amount  of  thrust  given  by  the  linear  actuator  at  the  leg;  and  the  body  attitude 
is  controlled  by  applying  torque  at  the  hip.  A  state  machine  permits  switching  between 
states.  Subsequently,  Raibert  extended  the  same  mechanism  and  controls  to  a  3D  version 
of  the  monopod,  Fig.  1.2b,  and  to  biped  and  quadruped  robots  with  prismatic  legs. 


(a)  (b) 


Figure  1.2:  Raibert  planar  hopper  (a)  and  3-dimensional  hopper  (b).  From  M.  Raib¬ 
ert,  Legged  Robots  That  Balance ,  MIT  Press,  1986 

The  Raibert  hopper  opened  the  door  to  a  long  series  of  SLIP-based  robots.  Zeglin’s 
Bow  Leg  hopper  [11]  is  composed  of  a  leg  controlled  with  strings  that  store  energy 
during  flight  but  are  loose  during  ground  contact,  so  that  it  behaves  like  a  spring- mass 
system.  A  planar  hopper  mounted  on  a  boom,  the  Bow  Leg  hopper,  has  been  also 
extended  to  a  3D  version  ([36]  and  [37]).  An  autonomous  and  untethered  hexapod  robot, 
RHex,  is  proposed  by  Saranli  et  al.  [38].  With  its  6  compliant  legs,  this  robot  is  able  to 
undergo  uneven  terrains  with  obstacles  exceeding  its  body  clearance,  with  limited  energy 
consumption.  The  ARL-Monopod  II  [26]  (and  a  previous  version  [39])  is  a  hopping  robot 
with  a  prismatic  leg  attached  to  the  body  via  a  torsional  joint  at  the  hip,  hopping  on  a 
treadmill.  Jonathan  Hurst’s  BiMASC  [40]  and  the  most  recent  ATRIAS  robot  have  been 
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mechanically  designed  to  match  the  dynamics  of  the  active  SLIP  model,  with  particular 
focus  on  energy  efficiency  (see,  for  example  [41]). 

1.2  Motivation  and  Contributions 

With  so  much  pre-existing  research  on  the  SLIP  model,  one  can  ask  “What  is  still 
left  to  be  done?” .  Surprisingly,  while  most  SLIP-like  robots  are  equipped  with  compliant 
legs,  there  is  little  investigation  on  the  effect  that  actuation  has  on  the  system’s  dynamics. 
As  shown  in  the  previous  section,  actuation  is  mostly  applied  in  a  “blind”  way  (i.e. ,  via 
feedback)  or  by  mimicking  animals’  muscles.  Besides  inputting  energy  into  the  system 
to  compensate  from  losses  at  impact  or  due  to  friction,  actuation  is  employed  to  reach 
a  certain  state  in  terms  of  running  speed/height  or  stride  length.  What  is  missing  is  an 
understanding  of  how  a  different  amount  of  actuation  applied  at  different  times  during 
stride  can  modify  the  state  reached  by  the  robot.  In  particular,  is  it  possible  to  use 
actuation  to  simplify  the  system’s  dynamics?  What  is  the  reachable  state-space  of  the 
system,  and  what  actuation  policy  enables  the  system  to  reach  the  most  states? 

The  goal  of  this  thesis  is  to  explore  how  a  series  elastic  actuator  applied  to  the 
SLIP  model’s  leg  can  change  the  system’s  dynamics.  We  principally  investigate  how, 
through  actuation,  we  can  solve  partially  or  completely  the  system’s  equations  of  motion, 
and  enforce  a  desired  trajectory.  We  also  determine  the  reachable  state-space  of  the 
model  using  several  different  actuation  strategies,  and  we  investigate  the  variation  of  the 
reachable  state  with  respect  to  the  particular  actuator  motion,  providing  a  relationship 
between  direction  of  actuator  displacement  and  location  of  the  reached  state.  Since  the 
SLIP  model  can  be  seen  as  a  template  for  more  complex  systems,  it  is  realistic  to  believe 
that  studies  on  the  effect  of  actuation  can  be  applicable  to  more  realistic  dynamics. 

This  thesis  is  structured  as  follows.  Chapter  2  introduces  the  passive  and  active  SLIP 
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models,  characterizes  their  dynamics,  underlines  the  differences  between  the  two,  and 
contains  symbols  and  definitions  that  will  be  used  throughout  the  document.  Chapter  3 
investigates  the  shape  of  the  reachable  space  for  the  active  SLIP  model  for  different 
actuator  motions.  Chapter  4  proposes  an  actuation  strategy  to  enforce  a  symmetric 
trajectory  for  the  stance  phase,  and,  at  the  same  time,  to  allow  us  to  compute  an  ana¬ 
lytical  solution  for  the  system’s  dynamics.  In  Chapter  5  we  propose  a  partial  feedback 
linearization  action  for  actuator  displacement  to  analytically  solve  part  of  its  dynamics, 
thereby  reducing  computational  time  and  increasing  the  practicality  of  performing  on¬ 
line  control  actions.  This  is  then  paired  with  a  two-part  control  action  to  add/remove 
energy  to/from  the  system  and  modify  the  upcoming  apex  state  to  span  an  open  set 
within  the  reachable  apex  states.  In  addition,  we  develop  two  control  strategies  for  on¬ 
line  computation  of  actuator  displacement  and  leg  positioning:  one  to  drive  the  system 
to  a  desired  state,  even  in  the  presence  of  terrain  perturbation;  the  other  to  control  the 
system  to  hop  on  a  desired  set  of  terrain  footholds.  Furthermore,  we  propose  an  adap¬ 
tive  control  technique  for  steady-state  locomotion  on  flat  terrain  to  reduce  computation 
errors  by  the  use  of  an  approximation  of  the  leg-angle  dynamics  during  the  stance  phase, 
and  we  demonstrate  the  proposed  strategy  on  a  more  dynamically  sophisticated  planar 
hopper  model.  In  Chapter  6  we  investigate  how  different  actuator  motions  can  affect 
the  system’s  state,  and  we  propose  a  control  strategy,  based  on  graphical  and  numerical 
studies  of  the  reachability  space,  and  updated  throughout  the  stance  phase,  to  drive  the 
system  to  a  desired  state.  We  quantify  its  performance  benefits,  particularly  in  serving 
as  an  error-recovering  method.  The  objective  of  our  control  strategy  is  not  to  replace 
any  leg-placement  approach  proposed  by  other  works,  but  rather  to  be  paired  with  any 
other  leg-placement  or  path  planning  method.  Its  main  advantage  is  the  ability  to  reduce 
the  effects  of  sensing  errors  and  disturbances  happening  at  landing  as  well  as  during  the 
stance  phase.  Finally,  Chapter  7  contains  conclusions  and  future  work. 
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The  SLIP  Model:  Preliminaries 


In  this  chapter  we  review  the  structure  of  the  passive  SLIP  model  and  its  dynamics.  We 
also  introduce  its  actuated  version  that  will  be  used  throughout  the  thesis,  the  so-called 
active  SLIP  model,  explaining  in  detail  its  dynamics  and  the  modifications  made  to  the 
original  model  to  incorporate  energy  variations. 

2.1  Passive  SLIP  Model 

The  passive  SLIP  is  modelled  as  a  point  mass,  M,  attached  to  a  massless  spring  leg, 
with  length  l  and  spring  stiffness  constant  k ,  as  shown  in  Fig.  2.1.  The  terminal  part 
of  the  leg  is  referred  to  as  the  foot.  Running  dynamics  for  the  SLIP  model  occur  on 
the  sagittal  plane,  and  constitute  an  hybrid  system.  In  fact,  they  consist  of  two  phases 
(see  Fig.  2.2a):  the  flight  phase ,  where  the  body  is  in  the  air  and  follows  a  ballistic 
trajectory;  and  the  stance  phase,  where  the  foot  is  in  contact  with  the  ground,  and  the 
compression/extension  of  the  spring  completely  defines  the  mass  dynamics.  We  will  call 
touch-down  (TD)  the  instant  that  marks  the  transition  from  flight  to  stance;  and  take-off 
(TO)  the  instant  that  marks  the  transition  between  stance  and  flight.  During  the  flight 
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phase  gravity  is  the  only  force  acting  on  the  system.  Defining  as  x  and  y  respectively  the 
forward  and  vertical  coordinates  of  the  mass,  the  equations  of  motion  during  flight  can 
be  written  as: 


x{t)  =  0, 
y(t)  =  -g, 

where,  as  customary,  g  is  the  gravitational  acceleration.  The  highest  point  reached  by 
the  mass  during  flight  is  called  the  apex  state ,  and  it  is  defined  by  zero  vertical  velocity, 
i.e. ,  ya  =  0  [m/s].  Therefore,  the  apex  state  is  completely  defined  by  a  three-dimensional 
vector  s  =  {xa,  ya ,  4}- 

The  stance  phase  starts  with  the  leg  hitting  the  ground  with  a  touch-down  angle 
9td •  While  the  body  moves  forward,  the  spring  compresses  until  it  reaches  its  minimum 
compression  point,  then  it  starts  expanding  (see  Fig.  2.2b).  When  the  spring  reaches 
its  equilibrium  position  (i.e.,  when  the  forces  in  the  spring  are  back  to  zero),  the  system 
leaves  the  ground  with  a  certain  take-off  angle  Qto-  As  shown  in  Fig.  2.1,  we  define  £{t)  as 
the  leg-length  as  a  function  of  time,  and  9{t)  as  the  leg-angle  measured  counterclockwise 
with  respect  to  the  positive  horizontal  axis,  while  4  is  the  spring  length,  and  4,o  is 
the  spring  length  at  equilibrium.  The  state  of  the  mass  can  be  easily  converted  from 
Cartesian  coordinates  into  polar  coordinates: 

x(t)  =  £{t)  cos  9(t),  y(t)  =  £{t)  sin  9(t). 

Then,  the  Lagrangian  for  the  system  during  stance  phase  is  derived  as  follows: 

M  k 

L  =  t2e 2  +  £2)  -  Mg£ sin 9  -  -(4  -  4,o)2, 
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and  the  equations  of  motion  in  polar  coordinates  for  the  stance  phase  can  be  written  as: 

t  =  -t?(4  -  4,o)  -  gsinO  +  £02,  (2.1) 

IVl 

0  =  -2^0-|cos0.  (2.2) 


Note  that,  despite  their  simplicity,  equations  (2.1)-(2.2)  are  not  analytically  solvable. 

Let  us  introduce  the  commonly  used  non-dimensional  relative  spring  stiffness,  7,  that 
will  be  used  throughout  this  thesis,  defined  as: 


7  = 


(2.3) 


The  introduction  of  the  relative  spring  stiffness  is  motivated  by  the  fact  that,  if  we  focus 
on  biped  and  quadruped  robots  where  each  leg  functions  and  behaves  the  same  way  as 
the  other (s),  we  can  study  their  performance  by  considering  the  relative  spring  stiffness 
7  per  leg  (or  net  per  multiple  legs  in  contact). 

In  this  work,  we  will  call  jump  the  transition  from  one  apex  state  to  the  next.  A 
jump  is  successful  if  the  forward  velocity  during  the  entire  trajectory  only  takes  positive 
values.  Furthermore,  to  allow  the  leg  to  swing  forward  in  flight,  we  require  the  distance 
between  ya  and  the  terrain  height  to  be  bigger  than  the  leg  length  at  equilibrium  £0. 

\ 

\ 

Figure  2.1:  Classic  SLIP  model  (passive  SLIP) 


For  each  forward  velocity  (or  apex  height,  equivalently,  since  energy  is  conserved) 
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(b) 


Figure  2.2:  (a)  Scheme  that  illustrates  the  phases  of  the  SLIP  model  trajectory,  (b) 
Stance  phase:  the  point  of  maximum  compression  corresponds  to  l  =  0. 


there  is  a  unique  touch-down  angle  the  gives  a  symmetric  step,  i.e. ,  the  apex  state  is 
preserved  due  to  zero  net  acceleration  during  stance.  The  corresponding  point  on  the 


ground  is  called  the  neutral  point.  Variations  from  the  neutral  point  have  the  effect  of 


increasing  or  decreasing  the  next  forward  velocity.  For  forward  hopping,  if  the  foot  is 
positioned  before  the  neutral  point,  than  there  is  a  net  acceleration,  while  if  the  foot  is 
positioned  after  the  neutral  point  there  is  a  net  deceleration. 


2.2  Active  SLIP  model 


The  classic  SLIP  model  is  energetically  conservative:  the  repositioning  of  the  leg  to 
the  desired  touch-down  angle  does  not  require  any  energy,  and  no  energy  is  lost  during 
impact  with  the  ground  as  well.  However,  various  studies  on  legged  locomotion,  e.g.,  [8] 
and  [42],  suggest  that  legs  store  and  dissipate  energy  during  motion.  Based  on  this 
evidence,  we  modify  the  passive  SLIP  by  adding  to  the  leg  a  piston-like  actuator  in  series 
with  the  spring,  as  shown  in  Fig.  2.3.  We  will  refer  to  the  actuated  SLIP  model  as  the 
active  SLIP. 
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Figure  2.3:  Active  SLIP  model  (active  SLIP) 

Define  £act(f)  to  be  the  actuator  length,  and  dact, o  be  its  nominal  position  at  rest. 
Throughout  the  stance  phase,  the  actuator  can  continuously  extend  ( dact  >  £act, o)  °r 
retract  ( £act  <  £act, o)  front  its  nominal  position  within  a  certain  displacement  range. 
Then,  the  equation  that  describes  the  evolution  of  the  leg  length  (2.1)  can  be  re-written 
as: 

£(t)  =  -  4-4ct(£))  -  gsind  +  W2,  (2.4) 

where  £q  =  £act,o  +  £k,o •  The  actuator  manages  absorption  and  production  of  energy 
during  the  stance  phase  by  compressing  and  decompressing  the  spring,  with  the  main 
advantage  of  allowing  the  energy  at  the  beginning  and  at  the  end  of  the  stance  phase 
to  be  different.  Positive  and  negative  values  of  £act  correspond  to  a  compression  and 
extension  of  the  spring,  respectively. 
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Given  an  inital  apex  state  s0  =  {xai  ya,  xa},  let  us  define  the  following  function  as  the 
transition  from  one  apex  to  the  next: 

o,  9td)  Si,  (3.1) 

in  the  case  of  the  passive  SLIP,  or,  for  the  active  SLIP: 


Xact(s  0,  9tDi  (■act  (t)  )  —  Si,  (3.2) 

where  Si  is  the  apex  state  reached  in  one  jump  from  the  initial  apex  state  Sq,  with  touch¬ 
down  angle  6TD  and  actuator  value  function  during  stance  phase  ^act(t)  (if  applicable). 

Then,  we  call  the  reachable  space  of  s0  the  set  of  all  apex  states  that  can  be  reached 
from  s o  in  one  jump  by  positioning  the  leg  and,  in  the  case  of  the  active  SLIP,  by  moving 
the  actuator: 

"^•pass(so)  •  {si  |  3  9tD  ■  fhpass(So,  O^Di  f-a(t) )  Si},  (3.3) 
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and 


T^act(s o)  {S1  |  3  0TD,  ^aif)  '■  <^act(sO)  @Td)  ~  sl}, 


(3.4) 


where  6td  C  [7t/2,  7t],  and  i act{t )  is  a  continuous  function  bounded  by  the  physical 
limitations  of  the  actuator,  such  as  limited  travel  and  velocity  and  torque  limits. 

For  the  passive  SLIP  model,  the  reachable  space  is  a  line  in  the  {x,  y,  x}— space,  or, 
when  fixing  the  touch-down  angle,  a  point  in  the  {x,  y,  x}—  space.  Adding  the  series 
actuator  has  the  effect  of  extending  the  reachable  set  to  a  3-dimensional  surface,  whose 
shape  is  determined  by  the  ’’shape”  of  the  function  lact(f)  during  the  stance  phase,  and  is 
a  function  of  the  actuator’s  characteristics,  such  as  maximum  and  minimum  compression 
and  extension,  and  maximum  velocity  and  acceleration  allowed.  To  take  into  account 
the  limitations  of  a  real-life  actuator,  we  model  the  actuator  dynamics  as: 


tact(t)  =  Caet(fi)  +  4«(*i)(<  -  U)  +  k<*Jt  -  Uf , 


(3.5) 


where  kacc  is  the  acceleration  constant,  and  tj  is  the  time  at  which  movement  first  occurs. 
When  the  actuator  reaches  its  maximum  velocity  allowed,  vrnax ,  it  stops  accelerating 
and  proceeds  with  a  constant  velocity  motion.  Additionally,  the  maximum  actuator 
displacement  is  assumed  to  be  10%  of  the  leg  length  at  equilibrium,  £0,  and  the  actuator 
is  assumed  to  start  its  motion  with  zero  velocity  and  displacement.  By  fixing  the  touch¬ 
down  angle,  the  reachable  space  for  any  actuator  movement  consists  of  a  2-dimensional 
manifold,  as  shown  in  Fig.  3.1a  and  3.1b.  The  apex  states  {x,  y,  x}  reachable  in  one  step 
from  a  fixed  touch-down  angle  form  a  2-dimensional  manifold,  suggesting  that  the  three 
coordinates  are  not  independent  from  each  other.  In  Fig.  3.2,  the  complete  reachable 
space  for  different  values  of  the  actuator’s  acceleration  constant,  kacc  are  computed.  As 
one  might  expect,  the  higher  the  acceleration,  the  bigger  the  reachability  space.  In 
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(a) 


(b) 

Figure  3.1:  Subplots  (a)  and  (b)  show  two  different  views  of  the  3-dimensional  reach¬ 
able  space  in  one  step  for  a  fixed  touch-down  angle.  As  we  can  see,  the  set  is  a 
2-dimensional  manifold. 


these  particular  examples,  initial  conditions  were  chosen  to  achieve  a  passive  symmetric 
jump,  i.e.,  Xact(so,  0)  =  xpass(s0,  0td)s o  with  no  actuator  movement.  Note  that, 
due  to  the  lack  of  closed-form  solution  for  the  stance  dynamics  of  the  SLIP  model,  the 
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reachable  spaces  were  numerically  computed.  This  operation  generally  requires  significant 
computational  effort;  hence,  it  is  not  normally  feasible  to  compute  in  real  time  the 
reachable  set  for  any  particular  initial  condition. 


Figure  3.2:  Reachability  space  on  the  {y,  x}— space  for  a  passively  symmetric  jump 
with  actuator  moving  at  different  acceleration  values,  given  a  maximum  allowed  ve¬ 
locity  of  vmax  =  0.5  [m/s],  and  maximum  acceleration  constant  of  kacc  =  10  [m/s2]. 
The  ’’corner”  points  are  the  points  reached  with  maximum  actuator  movement  in 
the  positive  (upper  corner)  or  negative  (lower  corner)  direction.  The  black  x  mark 
represents  the  apex  state  passively  reachable  (i.e.,  without  any  actuation).  From  this 
graph  we  can  see  that  the  higher  the  maximum  acceleration  allowed,  the  bigger  the 
reachability  space  is,  in  terms  of  area.  Initial  conditions  chosen:  7  =  10,  initial  apex 
state  sq  =  {0,  1.3  [m],  4  [m/s]},  touch-down  angle  0td  =  122.14  [deg]. 


In  general,  for  varying  touch-down  angle,  the  reachable  space  is  a  set  of  2-dimensional 
manifolds,  one  for  every  touchdown  angle  employed,  forming  a  fan-shaped  3-dimensional 
volume,  as  shown  in  Fig.  3.3a  and  3.3a. 
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(b) 

Figure  3.3:  Subplots  (a)  and  (b)  show  two  different  views  of  the  3-dimensional  reach¬ 
able  space  in  one  step  for  varying  touch-down  angles.  The  dark  grey  manifold  corre¬ 
sponds  to  the  reachability  space  in  Fig.  3.1a  and  3.1b. 
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Enforced  Symmetry  of  the  Stance 
Phase 


Let  us  assume  we  want  to  control  the  system  dynamics  to  follow  a  certain  trajectory 
during  stance.  From  equation  (2.4),  the  actuator  displacement  necessary  to  maintain  a 
desired  feasible  trajectory  9des(t)  and  £des(t)  during  the  stance  phase  is  given  by: 

Lctit)  =  ^-(tdesit)  +  gsmddes(t)  -  £des(t)9des(t ))  +  £des(t)  -  4-  (4.1) 

For  the  control  action  to  succeed,  the  commanded  trajectory  needs  to  be  feasible:  the 
trajectory  we  want  to  enforce  needs  to  be  consistent  with  the  dynamics  of  the  system 
and  its  limitations.  Since  the  dynamics  of  9(t)  and  £{t)  are  strongly  coupled,  it  is  not 
possible  to  control  them  separately  by  applying  actuation  at  the  spring.  Furthermore, 
to  compute  (4.1),  the  desired  trajectory  needs  to  be  fully  known.  Given  we  can  observe 
the  initial  conditions  at  touch-down,  one  starting  point  in  planning  is  to  exploit  the  fact 
that  it  is  often  possible  to  reach  the  same  conditions  at  take-off.  When  this  is  the  case, 
the  straightforward  answer  is  to  find  a  symmetric  trajectory  which  gives  a  closed-form 
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solution  for  the  stance  phase. 

4.1  Enforced  symmetric  gait 

The  control  action  required  to  enforce  a  symmetric  gait  to  the  SLIP  model  is  not 
unique.  However,  not  all  possible  symmetric  laws  can  be  solved  in  closed  form.  In  what 
follows,  we  start  by  observing  that  the  passive  natural  dynamics  of  the  angular  velocity 
with  respect  to  the  angular  position  has  a  shape  that  is  similar  to  a  part  of  a  sine  wave 
(see  Fig.  4.1c).  Therefore,  a  logical  approach  is  to  enforce  the  following  dynamics  to  the 
angular  velocity: 


0(t)  =  Hsin6)(t)  +  c±, 


(4.2) 


where  the  parameters  A  and  ci  are  uniquely  determined  such  that  position  and  velocity 


of  the  mass  at  the  start  of  the  stance  phase  match  the  ones  at  the  touch-down  state. 
From  (2.2)  and  the  time  derivative  of  (4.2),  we  have: 


_  21td6td  +  9  cos  6td 
t tdOtd  cos  Otd 

Ci  =  Otd  ~  A  sin  Otd  , 


where  Otd,  Otd,  £ td ,  P-td  are  position  and  velocity  at  touch-down  of  leg  angle  and  leg 
length,  respectively.  We  can  then  solve  (2.1)  to  find  a  closed-form  solution  for  the  leg 
length  dynamics: 


(4.3) 


with 
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It  is  important  to  note  that  the  leg  length  £{t)  can  never  take  negative  values.  Since  the 
minimum  value  of  £{t)  corresponds  to  the  case  6{t)  =  7t/2,  then 


q  Co 

+  ,  .  ~  >  0. 


A(A  +  Cl)  y/A+~c[ 


For  the  leg  angle,  we  can  solve  (4.2)  to  find  a  closed-form  solution  for  the  angle  dynamics 
during  stance: 

_!  / A  +  C4tanh(c4(t  +  c^)/2)' 


9{t)  =  —2  tan  1 


c  1 


for  i  6  [0,  txo\,  with 


and 


C3  =  -4tanh--('4  +  ‘:itan(^D/2) 
C4  V  C4 


Cl  =  JA2  -  cl 


(4.4) 


(4.5) 


The  take-off  time,  txoi  is  given  by: 


2  ,  -1  f 

tro  =  —C3 - tanh 

c  4  V 


_i  (A  +  ci  tan(6lTo/2)  j 


c4 


where  9to  is  flie  take-off  angle.  In  order  for  the  above  equations  to  be  feasible,  we  require 
the  following  constraint: 


\A  +  Ci  tan  {Qtd/ 2)|  <  | \J A2  —  cf | , 

and  the  same  constraint  can  be  enforced  for  all  9(t),  obtaining  the  following: 

9(t)  e  [2 tan-1  — — — j,  2 tan-1  ^ — — —  ^].  (4.6) 
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Note  that  the  parameter  C4,  defined  in  (4.5),  can  be  a  real  or  imaginary  number.  All 
equations  from  (4.4)  to  (4.6)  still  hold. 


Figure  4.1:  In  these  figures,  the  blue  dotted  line  refers  to  the  passive  SLIP,  while  the 
solid  red  line  refers  to  the  active  SLIP.  Figs,  (a)  and  (b)  show  the  mass  trajectory  for 
passive  and  active  SLIP  with  respect  to  the  horizontal  position  and  with  respect  to 
time,  respectively.  Fig.  (c)  shows  the  evolution  of  a  =  it  — 9  and  a  during  stance  for  the 
passive  and  active  SLIP,  and  Fig.  (d)  represents  the  required  actuator  displacement 
A£act.  The  parameters  used  are  k/M  =  500[s~2],  =  l[m]  and  0  =  0.5[m].  The 

initial  apex  state  is  sa  =  {0  [m],  1.4  [m],  2  [m/s]},  the  touch-down  angle  is  6td  =  100 
[deg] ,  and  terrain  height  is  0  [m] . 


Plots  in  Fig.  4.1  show  an  example  of  the  mass  trajectory  and  the  angular  acceleration 


24 


Enforced  Symmetry  of  the  Stance  Phase 


Chapter  4 


of  passive  and  active  SLIP  for  the  same  initial  conditions,  and  the  required  actuator 
displacement  for  the  active  case. 

4.2  Control  actions 

4.2.1  Fixed  leg  length  at  take-off 

When  the  trajectory  is  symmetric,  we  have  that  6to  —  7T~  @td ,  and  the  apex  velocity 
components,  xa  and  ya,  and  the  mass  apex  height  ya  are  preserved  from  one  apex  state 
to  the  next.  The  only  control  action  available  is  the  touch-down  angle  Qtd,  and  its  choice 
affects  the  forward  position  xa  at  the  next  apex  state.  If  we  start  from  the  same  initial 
apex  state  and  we  use  the  same  set  of  possible  touch-down  angle  Otd ,  it  is  clear  that 
enforcing  a  symmetric  trajectory  changes  the  set  of  reachable  apex  states  with  respect  to 
the  set  reachable  by  the  passive  SLIP.  As  an  example,  Fig.  4.2  shows  the  set  of  reachable 
apex  states  for  the  active  and  passive  SLIP.  We  can  see  that,  despite  the  different  spatial 
direction  of  the  two  curves,  the  span  of  reachable  points  is  comparable  in  size.  This 
suggests  that,  even  if  the  reachable  apex  states  are  not  the  same,  enforcing  a  symmetric 
gate  does  not  generally  represent  a  limitation  with  respect  to  the  passive  SLIP,  but  it 
rather  changes  the  direction  of  the  reachable  space. 

4.2.2  Variable  leg  length  at  take-off 

There  are  some  cases  when  enforcing  a  symmetric  gait  can  pose  a  limit  on  the  capa¬ 
bility  of  successfully  undergoing  a  particular  terrain.  For  example,  on  terrains  with  tall 
obstacles,  or  with  stairs,  it  may  be  necessary  to  add  (or  remove)  energy  to  the  system 
to  achieve  a  successful  jump.  Therefore,  we  allow  the  leg  to  leave  the  ground  either 
before  or  after  the  state  which  mirrors  the  touch-down  state  has  been  achieved.  This 


25 


Enforced  Symmetry  of  the  Stance  Phase 


Chapter  4 


Figure  4.2:  The  dotted  blue  line  and  the  solid  red  line  represent  the  reachable  points 
for  passive  and  active  SLIP,  respectively,  with  touch-down  angle  6td  €  [10,  90]  [deg]. 
The  system  parameters  are  k/M  =  500[s~2],  £q  =  1  [nr],  initial  apex  state 
sa  =  {0  [m],  1.4  [nr],  2  [m/s]},  and  terrain  height  0  [nr].  Note  that  both  lines  lie  on  a 
surface  of  constant  energy. 


simply  requires  setting  £act  such  that  the  leg  spring  is  fully  unloaded,  i.e.,  such  that 
£act  =  £  —  £o  +  Pactfi-  It  is  then  possible  to  choose  the  leg  length  at  take-off,  or  equiv¬ 
alently  to  choose  the  take-off  angle  Q^o  to  be  different  than  7r  —  9td-  Hence,  the  mass 
will  leave  the  ground  at  a  different  take-off  state,  ( xto ,  Uto ,  xto,  Vto ),  and  will  reach  a 
different  new  apex  state  than  otherwise  reached.  In  particular,  the  velocity  components 
at  take-off,  xto  and  yxoi  will  be  different  than  the  velocity  components  at  touch-down, 
xtd  and  Vtd-  Fig.  4.3  shows  the  set  of  all  the  apex  states  reachable  from  an  initial 
nominal  apex  state.  The  coloring  scheme  reflects  the  direction  of  the  next  apex  in  state 
space  with  respect  to  controlling  the  touch-down  angle  6td ,  Fig.  4.3a,  or  the  leg  length 
at  take-off  £toi  Fig.  4.3b. 
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Figure  4.3:  These  plots  represent  the  set  of  apex  states  reachable  from  the  ini¬ 
tial  apex  state  sa  =  {0  [m],  1.4  [m],  2  [m/s]}  in  the  actuated  case  with  leg  length  at 
take-off  Ito  €E  [0.5,  1.2]  [nr]  and  Otd  €  [100,  117]  [deg].  The  system  parameters  are 
k/M  =  500 [s  2] ,  ^o  =  l  [nr],  and  the  terrain  height  is  0  [nr]. 

The  main  advantages  of  such  a  procedure  are  twofold.  First,  the  capability  of 
adding/subtracting  energy  to/from  the  system  allows  us  to  reach  apex  states  otherwise 
unreachable,  which  translates  into  the  ability  of  operating  in  a  wider  variety  of  terrains. 
Second,  while  the  choice  of  Otd  acts  as  a  control  action  that  happens  before  the  impact 
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with  the  terrain,  the  choice  of  6to  is  done  after  the  impact.  This  represents  a  big  advan¬ 
tage  in  case  of  noisy  measurements  (e.g.,  on  the  terrain,  on  the  position  of  the  leg,  etc.), 
since  it  is  now  possible  to  partially  act  against  the  effect  of  the  noise  during  the  actuated 
stance  phase. 

4.3  Applications 

This  section  contains  two  applications  of  the  proposed  symmetric  enforcement  strat¬ 
egy.  Note  that  the  applications  are  manifold,  and  here  below  we  just  show  two  useful 
examples.  In  the  first  application,  we  assume  the  terrain  is  divided  into  allowed  and  not 
allowed  foothold  sets,  and  we  want  to  determine  a  path  from  an  initial  foothold  to  an  end 
foothold.  In  the  second  application,  we  aim  to  reach  a  desired  constant  forward  velocity 
or  constant  apex  height  during  motion  on  unknown  rough  terrain.  We  will  assume  that 
the  actuator  can  be  extended/contracted  only  up  to  a  certain  percentage  of  its  equilib¬ 
rium  length,  and  define  with  £™rfx  the  maximum  feasible  actuator  length.  This  will  pose 
a  constraint  on  the  velocity/height  reachable  in  one  step. 

4.3.1  Foothold  placement 

One  of  the  main  applications  of  the  closed-form  solution  is  the  ability  to  plan  a 
trajectory  given  a  set  of  feasible  footholds:  for  example,  trajectory  planning  on  a  terrain 
with  forbidden  areas,  such  as  holes  or  water. 

For  any  initial  apex  state  {xa,  ya,  xa}  and  terrain  height  h,  we  can  identify  a  minimum 
and  a  maximum  touch-down  angle,  9mm  and  9max ,  that  allow  the  system  to  complete  a 
successful  jump.  Let  us  call  step  length  the  distance  between  two  subsequent  landing 
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points  (Fig.  4.4),  and  let  us  define  the  two  quantities  dm  and  d^j  as 


dm  j' a  \  /  (]Ja  h  hj  sill  drn(lx)  f(j  COS  9 max  > 

9 


(4.7) 


Am  =  iot/  ~(Va  -  h  -  la  sin  9min)  -  t„  cos  0„ 


When  the  gait  is  symmetric,  the  minimum  and  maximum  step  lengths  at  each  apex  state 
{xa,  ya,  xa}  depend  on  the  length  of  the  previous  step  as  follows:  the  minimum  step 
length  is  sm  =  d0  +  dm,  and  the  maximum  is  %  =  d0  +  dM,  where  d0  is  the  distance 
between  xa  and  the  previous  landing  point.  It  follows  that  the  distance  between  two 


Figure  4.4:  The  figure  shows  the  distance  between  touch-down  locations  for  two  con¬ 
secutive  steps. 


feasible  footholds  has  to  take  values  between  [rlO  +  m,dO  +  M\.  It  is  straightforward 
to  show  that  do  is  bounded:  do  G  [dm,  o?m]-  Therefore,  a  necessary  but  not  sufficient 
condition  for  foothold  placement  is  that  the  distance  between  two  feasible  footholds  has 
to  be  between  [2dm,  2 d^-].  It  is  clear  that  the  required  distance  between  footholds  is  not 
an  absolute  measure,  but  it  depends  on  the  previous  step  length  and,  in  the  case  of  rough 
terrain,  on  the  terrain  height.  In  fact,  values  (4.7)  depend  on  the  terrain  height  at  the 
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landing  point,  h.  For  simplicity,  we  will  assume  the  terrain  to  be  flat;  however,  what 
follows  can  be  easily  adjusted  in  case  of  rough  terrain. 

Assume  now  that  the  terrain  is  divided  into  allowed  and  not  allowed  foothold  sets,  and 
that  we  want  to  determine  a  path  from  one  foothold  Vq  to  another  one  vend  (or  to  a  set  of 
footholds,  Vend),  given  an  initial  apex  state  and  initial  d0.  We  also  assume  that,  at  each 
step,  we  have  the  capability  to  see  the  terrain  up  to  a  certain  distance  D  from  our  current 
location.  Let  us  create  a  graph  where  each  vertex  represents  a  foothold.  Starting  from  an 
initial  foothold/vertex  vq,  and  initial  length  do,  if  D  >  do+m  we  add  an  edge  from  vq  to  all 
the  footholds  vt  whose  forward  distance  from  v0  is  d(v o,  v*)  G  [ d0  +  dm ,  min{.D,  d0  +  ^Af}], 
and  we  associate  to  the  edge  connecting  n0  and  vt  the  length  di  =  d{y 0,  Vi)  —  dQ.  (Note 
that,  if  D  <  do  +  dm,  there  is  no  feasible  foothold  for  that  specific  initial  condition.)  Then, 
for  each  vertex  vt,  we  repeat  the  same  process,  until  all  footholds  have  been  explored. 
As  output,  we  obtain  a  graph  for  each  feasible  path  from  node  Vo  to  node  v  G  Vend •  The 
optimal  path  can  then  be  chosen  based  on  the  specific  definition  of  cost  function:  e.g.,  we 
can  associate  a  cost  to  each  edge  based  on  the  length  covered,  or  based  on  the  probability 
that  the  foothold  location  is  correctly  sensed  to  be  feasible. 


Algorithm  1  compute  feasible  paths 
1:  V  =  {/o, ...  fn}  set  of  all  footholds 
2:  Vend  V  V  set  of  desired  path  end 
3:  Initial  state  v  =  /0,  initial  length  d0 
4:  graph  G  =  ( V ,  E,  W ),  with  V  =  {w},  E  =  0,  W  =  0. 

5:  for  Vi  G  V  s.t.  d(vi,  v )  G  min{D,  [d0  +  dm,  d0  +  dju]}  do 
6:  create  edge  ev^Vi  with  weight  di  =  d{v,  vt)  —  d0 

7:  cost  Ci 

8:  Id  :=  {V,  Vi},  Ei  :=  {E,  ev>Vi},  Hy  :=  {W,  d0} 

9:  Gi  :=  {Vi,  Ei  Wi} 

10:  Repeat  steps  5-10  for  v  =  Vi,  do  =  di,  and  G  =  Gi  until  no  more  possible  forward 

jump  is  feasible 

11:  end  for 

12:  return  all  graphs  G  that  end  in  Vend 
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The  complexity  of  Algorithm  1  grows  with  the  amount  of  knowledge  of  the  terrain 
available  and  the  number  of  feasible  footholds.  However,  we  can  intuitively  suspect  that, 
when  the  feasible  regions  of  the  terrain  are  high  in  number  and  equally  spaced,  a  complete 
knowledge  of  the  terrain  is  not  required.  Such  an  algorithm  is  practical  if  jump  transitions 
can  be  computed  quickly,  therefore  having  an  analytic  solution  is  a  necessary  condition 
to  keep  the  computation  complexity  low. 

4.3.2  Constant  forward  velocity/constant  apex  height 

Here,  we  consider  two  applications  for  our  variable  take-off  leg  length  strategy,  con¬ 
cerning  the  problems  of  reaching  a  desired  running  speed  and  of  reaching  a  desired  apex 
height  with  respect  to  the  terrain.  In  particular,  the  problem  of  reaching  a  desired  run¬ 
ning  speed  has  already  been  studied  in  the  literature  with  different  assumptions  on  the 
model  (see  [25]  and  [43]). 

First,  let  us  consider  the  problem  of  reaching  a  desired  running  speed  at  the  next  apex 
state,  Xd,  on  rough  terrain.  Given  the  initial  touch-down  conditions  (9td  and  ?td),  we 
want  to  compute  the  take  off  angle  that  allows  us  to  obtain  the  desired  forward  velocity, 
i.e.,  we  need  to  solve  the  following  equation: 


Xd  —  &to  cos  9to  ~  Oto^to  sin  6to ■  (4-8) 

Define  the  variable  z  =  \J A  sin  9to  +  D .  After  some  tedious  calculations,  one  can  show 
that  solving  equation  (4.8)  is  equivalent  to  solving  for  z  the  following  equation: 

Ac2^5  +  (2cig  A  2 A2Xd)z2  +  (A3C2  —  Ac\c2')z  +  (2 A^ g  —  2c^g)  =  0, 

which  cannot  be  solved  analytically  but  can  be  easily  solved  numerically.  In  order  to 
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obtain  xd,  the  leg  needs  to  leave  the  ground  at  an  angle  0TO  =  vr  —  sin 


corre¬ 


sponding  to  a  certain  take-off  leg  length. 


Since  the  actuator  ability  to  contract  /extend  is  limited  to  £^r,  then  lmax  =  +£k- 


From  equation  (4.3)  we  obtain 


q  Co 

-| - J—  < 


A9to  \/~Q 


—  ^ maxi 


TO 


and  therefore 


By  defining 


z  > 


c2  +  \/c|  +  4  gi  max  M 


9  f 


(4.9) 


,  c2  +  \Jc\  +  4:g£max/A  z 2  -  ci 

^  - ,  c5  = 


’  ”,J  A  ’ 

max  y  1 

we  can  compute  the  maximum  forward  velocity  obtainable  in  one  step  as: 


Tmacr. 


‘  Az2 


+  £i)(_4  +  d£l 

z A  2  2 


(4.10) 


-i2c6)  +  ^(c^  +  l). 


The  minimum  forward  velocity  can  be  easily  computed  as  the  value  the  velocity  takes  at 
the  time  t  when  6{t)  =  7t/2. 

Now,  let  us  consider  the  problem  of  maintaining  a  constant  apex  height  with  respect 
to  the  ground.  This  application  is  particularly  useful  in  situations  where  the  terrain 
height  is  progressively  increasing  and/or  decreasing.  In  fact,  it  is  crucial  to  maintain  a 
certain  minimum  distance  from  the  terrain  to  avoid  collision  of  the  leg  with  the  ground. 
At  the  same  time,  it  is  important  to  keep  the  relative  apex  height  limited  to  avoid 
gaining  excessive  velocity  during  the  ballistic  descent,  to  avoid  “bottoming  out”  the 
passive  spring,  with  the  resulting  impossibility  to  perform  a  successful  jump. 
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Let  us  assume  no  prior  knowledge  of  the  terrain,  in  which  case  our  strategy  focuses  on 
keeping  the  height  of  the  next  apex  constant  with  respect  to  the  ground  at  the  previous 
landing  point,  which  is  known  a  posteriori.  Given  the  initial  touch-down  conditions  ( 9td 
and  fra)  and  the  terrain  height  at  the  previous  foothold,  yt,  we  want  to  compute  the 
take-off  angle  that  results  in  obtaining  the  desired  apex  height  A  yd  with  respect  to  yt. 
Note  that  A yd  —  ya  —  yt ,  where  ya  is  the  next  apex  height.  Solving  this  problem  is 
equivalent  to  solving  for  Q^o  the  following  equation: 

1  /.  .  \2 

A y<i  =  £to  sin  6 to  +  ^  to  sin  6to  +  ^to^to  cos  #toJ  •  (4-11) 

As  for  the  constant  velocity  case,  this  is  equivalent  to  solve  for  z  =  a/ A  sin  dTO  +  C\  a 
10th-order  polynomial  equation  in  z  (not  shown  here  because  of  its  excessive  length). 
Because  of  limitation  on  the  actuator  length,  we  obtain  again  constraint  (4.9)  on  z. 
Correspondingly,  we  can  compute  the  maximum  and  minimum  apex  height  obtainable 
in  one  step. 

4.4  Simulations 

In  what  follows  we  illustrate  through  simulations  the  performance  of  our  example 
control  strategies  as  discussed  in  4.3.1  and  4.3.2.  The  system  parameters  used  are  defined 
in  the  following  table: 


k/M  [s-2; 

9 

111  /  s 2 

^act,  0 

m] 

0 

m] 

100 

9.81 

0.5 

0.5 

4.4.1  Foothold  placement 

The  following  simulations  refer  to  4.3.1.  Algorithm  1  has  been  tested  on  a  rough 

terrain  of  length  50  [m] .  Knowledge  of  the  terrain  ahead  has  been  assumed  to  be  D  =  10 

33 


Enforced  Symmetry  of  the  Stance  Phase 


Chapter  4 


[m] .  The  terrain  height  has  been  built  as  a  normal  distribution  with  0  mean  and  standard 
deviation  a  =  .1.  The  terrain  step  length  has  been  chosen  to  be  0.6  [m].  The  apex  state 
is  sa  =  {1.3  [m],  1.4  [nr],  3  [m/s]}  and  d0  =  1.3  [nr].  In  order  to  use  the  algorithm,  we 
discretized  the  terrain  in  intervals  of  size  0.1  [m],  and  we  considered  to  be  unfeasible  all 
the  points  xt  with  height  yt(xt)  such  that  \yt(xt)\  >  0.2  [m],  and  all  the  parts  close  to 
a  change  of  height,  in  order  to  avoid  leg  collision  with  the  terrain  during  stance.  The 
cost  function  has  been  chosen  to  be  the  final  covered  distance.  The  output  is  shown  in 
Fig.  4.5.  The  green  line  and  red  line  represent  respectively  the  feasible  and  unfeasible 
parts  of  the  terrain,  the  blue  dotted  line  represents  the  mass  trajectory,  and  the  black  x 
marks  are  the  landing  points. 
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x  [m] 


Figure  4.5:  The  figure  shows  the  path  computed  by  algorithm  1  where  the  cost  function 
is  defined  to  be  the  final  covered  distance.  The  green  and  red  line  represent  respectively 
the  feasible  and  unfeasible  parts  of  the  terrain,  the  blue  line  represents  the  mass 
trajectory,  and  the  black  x  marks  are  the  landing  points. 
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4.4.2  Constant  forward  velocity/constant  apex  height 

We  now  demonstrate  through  simulations  the  performance  typical  of  the  control 
strategies  introduced  in  4.3.2. 

First,  let  us  examine  the  problem  of  reaching  a  desired  forward  speed.  Consider  the 
active  SLIP  model,  and  assume  the  actuator  on  the  leg  can  extend/contract  up  to  20%  of 
the  leg  equilibrium  length.  Let  us  assume  the  active  SLIP  travels  on  a  rough  terrain,  but 
that  it  is  not  provided  with  any  knowledge  of  it.  We  will  then  treat  the  terrain  as  if  it  were 
completely  flat.  We  also  fix  the  touch-down  angle  at  9td  =  100  [deg].  During  stance, 
we  compute  the  take-off  angle  to  reach  the  desired  velocity  by  solving  (4.8)  as  follows. 
If  the  maximum  forward  velocity  obtainable  in  one  step,  imax  computed  with  (4.10), 
is  xmax  >  id,  then  the  take-off  angle  9to  will  be  computed  using  (4.8).  If  imax  <  id, 
then  9to  will  be  chosen  to  reach  at  the  next  apex  forward  velocity  imax ,  and  therefore 
the  number  of  steps  necessary  to  reach  the  desired  steady-state  will  be  greater  than  one. 
Fig.  4.6a  shows  a  numerical  example,  with  initial  apex  state  sa  =  {0  [m],  1.3  [m],  2  [m/s]} 
and  desired  forward  velocity  id  =  2.5  [m/s].  As  we  can  see,  the  system  reaches  the  desired 
forward  velocity  in  2  steps,  and  it  is  able  to  successfully  run  on  a  rough  terrain  with  error 
that  often  exceeds  50%  of  the  leg  length.  Fig.  4.6b  shows  a  numerical  example,  with  initial 
apex  state  sa  =  (0  [nr],  2.5  [m],  2.3  [m/s]}  and  desired  forward  velocity  xd  =  2 [m/s], 
smaller  than  the  inital  apex  velocity.  As  we  can  see,  the  system  reaches  the  desired 
forward  velocity  in  1  step. 

Let  us  now  examine  the  problem  of  reaching  a  desired  height  from  the  terrain,  as¬ 
suming  the  actuator  can  extend/contract  up  to  20%  of  the  leg  equilibrium  length.  Let 
us  assume  the  active  SLIP  is  not  provided  with  any  knowledge  of  the  terrain  it  trav¬ 
els  on.  We  fix  the  touch-down  angle  at  9td  =  100  [deg].  During  stance,  we  compute 
the  take-off  angle  to  reach  the  desired  terrain  height  with  respect  to  the  previous  land- 
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(b) 


Figure  4.6:  The  figure  shows  two  numerical  examples  of  the  forward  velocity  control 
strategy  proposed  in  4.3.2.  The  black  line  represents  the  real  terrain,  and  the  dotted 
green  flat  line  the  sensed  terrain.  The  red  x  marks  are  the  landing  points  and  the 
dotted  blue  line  is  the  mass  trajectory. 
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ing  height  by  solving  (4.11).  Fig.  4.7a  shows  a  numerical  example,  with  initial  apex 
state  sa  =  {0  [m],  2.3  [m],  2.3  [m/s]}  and  desired  relative  height  A yd  =  1.5  [m],  whereas 
Fig.  4.7b  shows  a  numerical  example,  with  initial  apex  state  sa  =  {0  [m],  1.3  [m],  2  [m/s]} 
and  desired  relative  height  A  yd  =  1.7  [m]. 


(a) 


Figure  4.7:  The  figure  shows  two  numerical  examples  of  the  relative  height  control 
strategy  proposed  in  4.3.2.  The  black  line  represents  the  real  terrain,  the  dotted  blue 
line  is  the  mass  trajectory,  the  red  x  marks  are  the  landing  points  and  the  magenta 
dots  are  the  apex  states. 
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4.5  Conclusions 

In  this  chapter,  we  developed  a  control  law  for  actuator  displacement  in  order  to 
enforce  a  symmetric  stance  phase  to  the  model.  For  the  specific  control  law  we  enforced, 
we  were  able  to  compute  a  closed-form  solution  for  the  trajectory  of  the  mass  during 
the  stance  phase.  We  then  exploited  the  possibility  of  the  leg  leaving  the  ground  either 
before  or  after  reaching  the  take-off  state  that  is  symmetric  to  the  touch-down  state, 
and  as  a  consequence  we  were  able  to  add  and  remove  energy  from  the  system.  We 
also  provided  two  examples  of  useful  applications  for  our  method  (trajectory  planning 
on  fixed  footholds,  and  running  on  rough  terrain  at  constant  forward  velocity  and  at  a 
constant  relative  height),  and  we  validated  them  through  simulations. 
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Approximation  and  two-element 
control 


We  have  seen  in  Chapter  4  that  the  series  elastic  actuator  can  be  used  to  enforce  a 
trajectory  for  the  stance  phase.  In  particular,  the  trajectory  was  chosen  to  be  symmetric 
to  be  able  to  analytically  solve  the  equations  of  motion.  There  are  undeniable  benefits  in 
having  access  to  a  closed-form  solution  for  the  system’s  dynamics,  both  on  a  control  and 
on  a  computational  point  of  view.  Unfortunately,  it  is  not  always  practical  to  enforce  a 
symmetric  stance  phase.  Especially  if  the  ability  of  locking  the  leg  at  take-off  is  precluded, 
it  is  not  possible  to  inject/remove  energy  to/from  the  system. 

In  this  chapter  we  investigate  how  it  is  possible  to  use  actuation  to  analytically 
solve  part  of  the  stance  phase  dynamics.  Subsequently,  we  pair  this  solution  with  an 
approximation  for  the  remaining  equation,  and  we  propose  a  control  strategy  to  drive 
the  system  to  a  desired  apex  state. 
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5.1  Approximating  the  stance  phase  dynamics  through 

PFL 

Here,  we  propose  a  strategy  for  the  actuator  displacement  4ct(f)  to  solve  part  of  the 
stance  dynamics,  and  we  propose  an  approximation  of  the  remaining  equation. 

5.1.1  Exact  solution  via  partial  feedback  linearization 

Let  us  write  the  total  actuator  displacement  as:  4ct(f)  =  V-niif)  +  -4(f),  with  to¬ 
tal  velocity  vact(t)  =  vni(t )  +  vc(t).  The  first  term,  4z(f),  performs  a  partial  feedback 
linearization  (PFL):  it  has  the  purpose  of  cancelling  the  nonlinear  terms  in  (2.4): 

4z(f)  =  sin 0(f)  -  £(t)6(t)2).  (5.1) 

Note  that,  at  touch-down,  this  becomes: 

£nl,TD  =  —  (g  sin  9td  —  £o@Td)i 

which  is  in  general  not  zero,  and  therefore  the  spring  may  need  to  be  pre-compressed  or 
pre-extended  during  flight.  Substituting  (5.1)  in  (2.4),  we  obtain: 

m  =  -L«<)  -  4  -  4W), 

where  4(f)  is  chosen  as  follows. 

We  drive  the  second  term,  4(f),  to  a  particular  value  4-  However,  since  the  actuator 
does  not  move  instantaneously,  we  assume  that,  after  cancelling  the  nonlinearity,  the 
actuator  moves  with  a  constant  velocity  vc  from  its  initial  position,  4(f*)  ,  until  it  reaches 
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the  desired  value  lc\ 

4 i(t)  +  4 (4)  +  vct,  if  | vct  +  4(4)1  <  |4 

4  ct(t)  =  < 

£ni(t )  +  4>  otherwise, 

as  illustrated  in  Fig.  5.1. 


Figure  5.1:  Example  of  the  total  actuator  value  £act{t )  during  stance  phase,  and  its 
two  components  £n i(t)  and  4(4 

Note:  in  general,  vni(t)  is  not  a  constant  value.  Since  vc  is  set  to  be  a  constant,  the  total 
actuator  velocity  required  vact(t)  is  a  time- varying  function. 

During  the  transition  time  required  to  reach  the  desired  actuator  value  £c,  the  equation 
that  describes  the  leg-length  dynamics  can  be  analytically  solved  as 


£(t)  —  r  +  vct  +  a  cos  (cut  +  /?),  (5.2) 
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where1 


/3  =  atan2  (  — 


a  = 


?(ti) 


LU 

4(4) 


cos  f3 


t(ti)  -4  -  4 (4)), 
r  =  4  +  4(4), 


4  is  the  initial  time,  and  cu  =  \Jk/M. 

Once  the  actuator  reaches  the  desired  final  value  4,  the  leg-length  dynamics  are 
described  by 

£(t)  =  r  +  a  cos  (tut  +  f3),  (5.3) 


where 


/3  =  atan2(-^^-,  44)  -  4  -  4), 

cu 

44)  -  4-4  o  ,  a 
a  = - 7^ — ,  r  —  4  +  4, 

cos(p) 


and  4  is  the  time  at  which  vctc  =  4- 


5.1.2  Approximation  of  the  leg  angle  dynamics 

The  equation  of  motion  for  the  angular  displacement  over  time,  0(f)  (2.2),  is  still 
not  analytically  solvable.  However,  we  notice  that  Equation  (5.3)  has  the  same  form  as 
the  approximation  of  the  leg-length  dynamics  in  [34],  with  two  main  differences:  (z)  [34] 
considers  the  passive  SLIP  model  only,  and  (ii)  the  equation  for  the  leg-length  dynamics 
provided  in  [34]  is  an  approximation,  while  (5.3)  and  (5.2)  are  exact  solutions.  We  then 
chose  to  follow  the  same  initial  steps  of  the  procedure  to  approximate  0(f)  proposed 
in  [34],  modifying  and  extending  the  results  to  adapt  them  to  our  actuated  case. 

■'Por  any  point  (x,  y)  in  the  xy-plane  minus  the  origin,  arctan2(y,  x)  is  defined  to  be  the  angle  between 
the  horizontal  positive  axis  and  the  point  (x,  y )  measured  counterclockwise. 
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When  the  actuator  reaches  the  final  desired  value  £c,  the  dynamics  of  the  leg  length 
during  stance  are  described  by  Equation  (5.3).  Let  us  define  a(t)  =  9(t )  —  tt/2,  and 
assume  that  both  a(t)  and  the  angular  span  Act  =  oitd  —  (Xto  ate  sufficiently  small  to  be 
simplified  as  sin(a)  ps  a.  Then,  equation  (2.2)  can  be  simplified  and  re-written  in  terms 
of  a,  becoming: 


(5.4) 


Toward  solving  for  a  corresponding  analytic  approximation  for  a{t ),  we  introduce  the 


variables  u{t)  and  p{t)  such  as  a{t)  =  p(t)u(t).  Substituting  a(t)  in  (5.4),  we  obtain 


Setting  the  coefficients  associated  to  ii  to  be  equal  to  zero  yields 


£p  +  £p  =  0, 


whose  solution  is  p  =  l/£.  Then,  equation  (5.4)  becomes: 


(5.5) 


When  the  actuator  reaches  the  final  desired  value  £c,  the  equation  that  describes  the 
dynamics  of  the  leg  length  during  stance  is  Equation  (5.3).  Then,  (5.5)  becomes 
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with  r  =  £(tc)  +  £c,  and  e  —  We  can  write 

z 

e  = - 

cos  p 


where 

£(tc)  £(ti ) 

7  =  - or  7  =  - 

4  +  4’  ~  4  +  4(4)  ’ 


Let  us  assume  small  leg  compression  and  small  actuation,  i.e.,  (£(£)  —  4— 4)/(4  +  4)  *C  1, 
then  we  have  that 


, ,  ,  m  , 

lim  e  =  - - - =r-  . 

2— ta(4  "h  4) 

which,  for  typical  values  of  w  and  leg  length  velocity,  has  magnitude  less  than  1.  Then, 
for  small  values  of  e,  we  have  that 


— -  «  1  -  e  cos  (wt  +  P)  +  e2  cos2  (wt  +  £) - , 

1  +  e  cos  (uit  +  p) 

and  thereby  obtain  the  Mathieu  equation 


u  —  u  (A2  —  eS  cos  (ut  +  f3))  , 


with 


5  =  co2  +  g/r ,  A  2  =  g/r. 


Equation  (5.6)  can  be  developed  as 


u(t)  =  u0(t )  +  euiit)  +  e2u2(t)  -| - . 


(5.6) 
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We  then  obtain 


a(t) 


1 

W) 


(' u0(t )  +  CMi(t)), 


where  u0(t)  and  Ui(t)  are  solutions  of 


d2  \  2 
if-u  °  ~ A 


(it2 


Ml  —  A2Wi 


0, 

— Suq  cost at  +  f3. 


We  want  now  to  modify  the  procedure  to  approximate  the  leg  angle  dynamics  pro¬ 
posed  in  [34]  and  performed  above,  to  adapt  it  to  (5.2).  During  the  transition  time 
required  to  reach  the  desired  value  £c,  the  equation  that  describes  the  leg  length  dynam¬ 
ics  can  be  analytically  solved  as  (5.2).  Then,  equation  (5.5)  becomes 

-c2o;2  cos  (cat  +  f3)  +  g 

u(t)  —  u(t ) - 7 - — - =  0. 

ci  +  C2  cos  (c at  +  (3)  +  vct 

Defining  0  =  cut  +  f3,  we  have  4w(t)  =  and  ii(t)  =  u>2u"(ip),  where  (•)'  represents 

the  derivative  with  respect  to  0.  Then,  the  previous  equation  becomes: 


m"(0)  -  m(0)  (  —  1  + 


1  +  C  +  C0 

1  +  K  COS  0  +  £0 


=  0, 


(5.7) 


with 


f  = 


£c2<a  ig 

«  =  - >  C  =  - • 


Vr 


UWr 


Cita  -  /3vc 

We  can  expand  the  fractional  term  in  (5.7)  around  cos 0  =  0,  obtaining 


1  + 


1  +  C  +  £0 


1  +  K  COS  0  +  £0 
A2  —  5(e  cos  0  —  s2  cos2  0  +  e3  cos3  0  —  •  •  • ), 
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where 


Then 


-  _  2(k  -  f )  r2 

c  —  - — ,  A 

2  +  7t£ 


2C 

2  +  ^’ 


6 


K  2(1  +  0 

0  -  0  +  2  +  7 rf  ■ 


a(t) 


where  u0(t)  and  Ui(t)  are  solutions  of 


u'q  —  A  2u0  =  0, 

u[  —  \2ui  =  —Suocos'ip. 
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5.2  Choice  of  the  actuator  constant  value  £c:  two- 
element  strategy 

We  now  propose  a  strategy  for  the  choice  of  the  actuator  constant  value  lc. 

Let  us  divide  the  stance  phase  in  two  parts,  separated  by  the  point  of  maximal  leg 
compression:  a  first  part,  where  £(t)  <  0,  and  a  second  part,  where  £(t)  >  0.  Our  main 
control  action  consists  in  choosing  two  constant  values  for  £c:  one  for  the  first  part,  £ci, 
and  one  for  the  second  part  of  the  stance  phase,  as  shown  in  Fig.  5.2. 


Figure  5.2:  Proposed  two-element  actuator  motion  during  stance  phase 

Dissimilarly  from  Raibert’s  famous  single- legged  hopper,  which  adjusts  the  energy 
via  actuation  once  at  mid-stance,  we  propose  to  adjust  the  actuator  motion  twice  during 
stance.  This  approach  allows  us  not  only  to  regulate  the  system’s  energy,  but,  more 
particularly,  its  two  components  of  forward  velocity  and  apex  height.  As  Fig.  5.2  il¬ 
lustrates,  actuation  during  the  entire  stance  phase  allows  the  system  to  reach  a  wider 
variety  of  apex  states  than  would  otherwise  be  reached  with  mid-stance  actuation  only. 
Specifically,  using  a  single  (1-dimension)  thrust  parameterization  will  clearly  map  to 

only  a  1-dimensional  set  of  apex  states,  while  the  full  reachable  set  of  apex  states  is  an 
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approximately  2-dimensional  surface  (  [44],  [15]). 

Furthermore,  the  decision  of  limiting  our  control  to  a  two-parameter  (4i,  fc2)  choice 
for  the  actuator  displacement  as  opposed  to  a  time  dependent  function  (as,  for  example, 
in  [13],  [44],  and  [15])  has  been  dictated  by  the  purpose  of  keeping  the  system  as  simple 
as  possible,  without  much  loss  on  performance.  Fig.  5.3a  and  5.3b  provide  an  example  of 
how,  by  setting  only  two  actuator  values,  it  is  possible  to  reach  in  one  jump  a  wide  range 
of  apex  states,  influencing  both  apex  height,  apex  velocity  and  apex  forward  position,  in 
all  reachable  directions. 


Figure  5.3:  Apex  states  { y ,  x}  (a),  and  { x ,  y,  x}  (b),  reachable  in  one  jump,  with 
M  =  10  [kg],  k  =  1962  [N/m],  4  =  1  [m].  Initial  apex:  xapex  =  0  [m],  yapex  =  1.4 
[nr],  Xapex  =  2  [nr/s],  Otd  =  103  [deg].  Maximum/minimum  actuator  length:  ±0.05 
[nr].  Actuator  moves  with  velocity  vc  =  0.5  [m/s].  Each  solid  blue  line  corresponds 
to  a  different  actuator  value  during  the  first  half  of  the  stance  phase,  4i>  while  each 
dotted  red  line  corresponds  to  a  different  actuator  value  during  the  second  half  of  the 
stance  phase,  42-  The  black  lines  correspond  to  the  case  of  either  4i  =  0  or  42  =  0. 


Because  the  reachable  {x,  y,  x}—  space  obtained  by  regulating  the  actuator  value 

(Fig.  5.3b)  takes  the  form  of  a  2-dimensional  surface,  increasing  the  dimension  of  the 

reachable  state-space  requires  that  we  add  another  degree  of  freedom.  This  can  be  done 
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by  varying  the  touch-down  angle  9td-  Then,  the  reachable  space  becomes  a  non-convex 
set  in  the  3-dimensional  space,  as  shown  in  Fig.  5.4. 


Figure  5.4:  Apex  states  {x,  y,  x}  reachable  by  our  control  strategy  in  one  jump  from 
the  initial  apex  state  {0, 1.4  [m] ,  2  [m/s]},  for  varying  values  of  9td  £  [85,  170]  [deg] 
(colorbar),  tc\  and  £.&  €  [—.05,  0.5].  Actuator  moves  with  velocity  vc  =  0.5  [m/s]. 
Relative  spring  stiffness  7  =  20,  and  leg  length  Iq  =  1  [m]. 


As  stated  earlier,  moving  the  actuator  results  in  a  change  of  the  system’s  energy.  The 
energy  at  apex  is  defined  as 

Ea  =  Mgya  +  ^Mx2a.  (5.8) 

Fig.  5.5  shows  the  change  in  energy  from  one  initial  apex  state  to  the  next  for  different 
values  for  the  couple  £&).  As  a  general  rule  of  thumb,  the  maximum  energy  increase 
is  obtained  by  extending  the  spring  in  the  first  half  of  the  stance  phase  (i.e. ,  ic\  <  0)  and 
compressing  it  during  the  second  half  (i.e.,  £&  >  0),  while  compressing  first,  and  then 
extending  (i.e.,  lc\  >  0  and  £&  <  0)  results  in  a  maximum  energy  decrease. 
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Figure  5.5:  The  contour  plot  represents  the  energy  change  values  A E  =  Ea 2  —  Ea\ . 
Ea) i  is  the  energy  (5.8)  at  the  initial  apex  state  ya  =  1.4  [m]  and  xa  =  2  [m/s].  Ea^ 
is  the  energy  (5.8)  at  each  next  apex  state  computed  for  Otd  =  103  [deg]  and  varying 
values  of  the  couple  (£c i,  lc 2). 


The  lack  of  body  inertia  in  this  simplified  model  prevents  the  use  of  body  attitude 
as  additional  control  action  to  reach  a  desired  state.  However,  it  is  not  unreasonable  to 
believe  that  the  strategy  we  propose  can  be  easily  paired  with  other  leg  placement  and 
body  attitude  control  strategies.  For  example,  the  leg  actuator  proportional  controller 
proposed  in  [12]  could  be  replaced  by  a  two-part  thrust  actuation  strategy  such  as  ours, 
preserving  the  energy  efficient  hip  actuator  controller  based  011  the  hip  passive  oscillations. 


5.3  Performance 


We  study  the  performance  of  our  monoped  hopper  in  terms  of  the  relative  spring 

stiffness  7  as  defined  in  (2.3).  Simulations  are  conducted  for  7  e  [10,  200],  using  constant 

values  for  i 0  and  M,  and  varying  k.  The  initial  apex  height  and  velocity  have  been 
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chosen  as  a  function  of  the  leg  length  and  the  time  scale  r  =  1  [s]  to  be  ya  G  [to,  2.5£0] 
and  xa  G  [0.5-^,  3^],  while  touch-down  angle  has  been  chosen  as  6td  G  [90,150]  [deg]. 
The  spring  length  at  equilibrium  has  been  assumed  to  be  £k$  =  0.5to,  with  a  maximum 
compression  of  £k,min  =  0.05£o- 

5.3.1  Feasibility 

It  is  important  to  point  out  that,  depending  on  the  system’s  parameters  and  its  initial 
conditions,  the  actuator  displacement  required  to  cancel  the  nonlinear  terms,  £ni{t),  could 
exceed  the  maximum  actuator  displacement  and  velocity  allowed,  or  could  bottom-out 
the  spring.  To  include  such  feasibility  constraints  in  our  work,  we  assume  the  total 
actuator  displacement  £act(t)  =  £c(t)  +  £ni{t)  must  not  exceed  10%  of  the  leg  length  £0, 
the  maximum  velocity  vact  =  vc  +  vni  not  exceed  £q /t,  and  that  £act{t )  at  any  given 
time  must  not  bottom-out  the  spring,  i.e.,  £k(t)  >  £k,min ■  The  allowable  amount  of 
displacement  and  velocity  for  the  nonlinear  part  {£ni)  and  the  constant  part  (£c)  can 
be  allocated  in  an  infinite  number  of  ways.  For  example,  Fig.  5.6a  shows  the  initial 
apex  states  that  require  £ni  <  0.5£o  and  vni  <  0.5£q/t  to  perform  a  symmetric  jump, 
while  Fig.  5.6b  considers  apex  states  that  require  £ni  <  £)A£0  and  vni  <  0.3£o /t.  The 
simulations  have  been  computed  for  several  values  of  7. 

Our  strategy  works  well  for  values  of  7  >  100  for  both  allocations  considered.  For 
7  <  10  the  size  of  the  set  of  feasible  initial  conditions  is  small,  posing  a  heavy  limit  to 
the  application  of  our  controlling  strategy. 

5.3.2  Error  reduction 

Now,  we  want  to  test  the  benefits  (in  terms  of  approximation  error)  of  our  proposed 
strategy  for  cancelling  the  nonlinear  terms  via  our  active  SLIP  control.  We  introduce  the 
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Figure  5.6:  Subplot  (a)  shows  the  initial  apex  states  that  require  lni  <  0.5f?o  and 
vni  <  0.5£o/t  to  perform  a  symmetric  jump,  while  subplot  (b)  considers  apex  states 
that  require  £ni  <  0.4f?o  and  vni  <  0.3£q/t.  Different  shades  refer  to  different  intervals 
of  7:  7  >  10  (black),  7  >  20,  7  >  40,  7  >  100,  and  7  >  200  (light  yellow). 


non-dimensional  percentage  errors  of  variables  x,  y  and  x,  respectively,  as: 


PEX  =  10011  X  B  X  l|2,  PEy  =  10011  V  B  V  Il2: 


PE±  =  10011  ^  a  X  2r, 


(5.9) 


where  y  and  x  are  height  and  velocity  at  apex  computed  via  approximation,  while  y  and 
x  are  the  actual  apex  height  and  velocity  computed  using  Matlab  numerical  solver  ode45, 
with  absolute  and  relative  tolerances  set  at  10~8.  The  time  constant  r  has  been  chosen 
to  be  equal  to  r  =  1  [s]. 

First  of  all,  why  is  it  useful  to  cancel  the  nonlinear  terms,  i.e. ,  what  is  the  benefit 
of  having  an  exact  solution  for  £{t)7  We  answer  this  question  by  comparing  the  per¬ 
centage  errors  (5.9)  for  the  approximation  proposed  in  [34]  versus  our  approximation 
with  nonlinearity  cancellation  (5.1).  Since  the  approximation  in  [34]  does  not  consider 
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actuation,  the  comparison  is  performed  with  respect  to  our  approximation  computed 
with  dact(t)  =  dni(t),  i.e.,  £c(t)  =  0.  Fig.  5.7a,  5.7b  and  5.7c  show  the  mean  percentage 
errors  PEX:  PEy  and  PE±  for  symmetric  and  non-symmetric  trajectories,  with  values 
of  7  G  [20,  200]  (values  of  7  <  20  have  not  been  considered  due  to  their  limitation,  as 
shown  in  Subsection  5.3.1).  We  can  see  that  our  proposed  strategy  significantly  reduces 
the  percentage  errors,  especially  for  lower  values  of  7.  This  can  serve  as  a  starting  point 
for  the  choice  of  7  while  building  a  hardware  prototype. 

We  now  compute  the  percentage  errors  (5.9)  for  the  actuated  SLIP  model,  with  our 
proposed  actuator  displacement  strategy  f?act(t)  as  in  (5.1.1).  Fig.  5.8a,  5.8b  and  5.8c  show 
the  mean  and  standard  deviation  of,  respectively,  the  percentage  errors  PEX,  PEy  and 
PE±,  computed  for  several  values  of  7  e  [20,  200]  and  a  set  of  60,000  initial  conditions. 
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Figure  5.7:  The  blue  stars  represent  the  PEs  computed  with  respect  to  the  stance 
phase  approximation  proposed  in  [34],  The  red  diamonds  represent  the  PEs  com¬ 
puted  using  our  proposed  approximation  via  nonlinearity  cancellation.  Initial  apex 
conditions  have  been  chosen  to  be  y  €  [£o,2.5£q],  x  €  [0.3^-,  3^],  and  9td  €  [85, 150] 
[deg].  Maximum  actuator  displacement  and  velocity  are  max(\£act\)  =  0.05^o  and 
max{\vact\)  =  0.5  £0/r.  ^ 
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Figure  5.8:  Percentage  errors  PEX  (a),  PEy  (b)  and  PE±  (c)  for  our  proposed  ap¬ 
proximation  with  actuator  displacement  (5.1.1).  The  colored  symbols  represent  the 
mean  values  of  the  percentage  errors,  and  the  vertical  bars  the  respective  standard 
deviations,  computed  for  a  pool  of  60,000  initial  conditions.  The  initial  apex  con¬ 
ditions  have  been  chosen  to  be  y  £  [4,  1.84],  x  €  [0.5^,  3^],  and  9td  €  [90,135] 
[deg].  Maximum  actuator  displacement  and  velocity  are  max(\£act\)  =  0.14)  and 

max(\vact\)  =  1  £q/t,  with  ic  €  [—0.054)  0.054],  and  vc  €  [0.54/t,  0.54/4-  If  either 

the  total  actuator  displacement  or  velocity  exceeds  the  maximum  values  allowed,  the 
actuator  is  assumed  to  saturate  during  numerical  computation  of  the  stance  phase 
trajectory. 
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5.3.3  Robustness  to  sensor  noise 

The  exact  solution  of  the  leg  length  dynamics  via  partial  feedback  linearization  and, 
more  generally,  the  approximation  of  the  overall  stance  dynamics,  depend  on  the  accu¬ 
racy  of  the  sensors’  measurements  available.  We  test  the  percentage  errors  (5.9)  for  the 
actuated  SLIP  model  when  sensor  noise  is  present.  In  particular,  we  consider  errors  in  the 
detection  of  the  initial  apex  height,  position  and  forward  velocity,  which  translates  into 
errors  of  the  stance  phase’s  initial  conditions,  and  noise  on  the  partial  feedback  lineariza¬ 
tion  term  (5.1).  Fig.  5.9a,  5.9b,  and  5.9c  show  the  mean  percentage  errors  when  sensor 
noise  is  present,  and  compare  it  with  the  noiseless  case.  Sensor  noise  on  the  initial  apex 
state  is  taken  from  a  uniform  distribution  within  ±1%  the  actual  values.  Furthermore, 
we  assume  the  controller  is  updated  at  1000  Hz,  and  the  partial  feedback  linearization 
term  (5.1)  is  affected  by  measurement  noise  modelled  as  a  Gaussian  distribution  with 
zero  mean  and  standard  deviations  a  —  0.1  for  0  and  6,  and  a  =  0.01  for  t.  While  the 
presence  of  noise  degrades  the  accuracy  of  the  approximation  of  the  system’s  dynamics, 
the  noise  considered  still  results  in  a  good  level  of  accuracy. 
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Figure  5.9:  Percentage  errors  PEX  (a),  PEy  (b)  and  PEX  (c)  for  our  proposed  approx¬ 
imation  with  actuator  displacement  (5.1.1),  with  and  without  sensor  noise.  The  blue 
diamond  dotted  lines  represent  the  mean  values  of  the  percentage  errors  with  noisy 
sensors,  while  the  red  star  solid  lines  represent  the  mean  values  of  the  percentage  errors 
with  perfect  sensing,  for  varying  values  of  7.  PEs  are  computed  for  a  pool  of  60,000 
initial  conditions,  chosen  to  be  y  €  [4,1-84],  *  €  [0.5^-,  3^r],  and  Otd  €  [90,135] 
[deg].  Maximum  actuator  displacement  and  velocity  are  max(\£act\)  =  0.14)  and 
max(\vact\)  =  14/r,  with  tc  €  [—0.054,  0.054],  and  vc  €  [0.54 A",  0.54 A"]-  Sensor 
noise  is  modelled  as  follows.  Noise  on  the  initial  apex  states  are  taken  from  a  uniform 
distribution  within  ±1%  of  the  actual  values,  while  the  partial  feedback  linearization 
term  (5.1)  noise  is  modelled  as  a  Gaussian  distribution  with  zero  mean  and  standard 
deviations  a  =  0.1  for  6  and  9,  and  a  =  0.01  for  £. 
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5.3.4  Cost  of  Transport 

At  each  time  t  during  the  stance  phase,  it  is  possible  to  compute  the  total  work  done 
by  the  spring  as: 

Wspring(t)  =  f  k(£k(t)  -  £k,0)i(t)dt.  (5.10) 

■'  tTD 

Note  that  £{t)  =  £k(t)  +  £act(t),  where  £act(t)  is  the  actuator’s  velocity. 

The  energy  stored  in  the  spring  is 

Ek{t)  =  —  (4(4  —  4,o)2- 

If  there  is  no  actuation,  the  energy  stored  in  the  spring  is  equivalent  to  the  work  done  by 
the  spring,  i.e. ,  Ek(t)  =  Wspring (t).  However,  when  actuation  is  present,  it  is  responsible 
for  part  of  the  energy  stored  in  the  spring:  Ek(t)  =  Wspring(t)  +  Wact(t),  where  Wact(t)  is 
the  work  done  by  the  actuator.  Then 


Wact(t)  =  Ek(t )  -  Wspring(t ) 


=  |(4(t)-4,0)2-  T  k(ek(t)  -  ekfi)e(t)dt 

^  J  tTD 


' tTD 


Pact{t)dt  +  — (4(4d)  —  4,o)2; 


and  Pact, 4 )  —  ~k(£k{t)  ~  4,o)4ct(4  is  ^ie  power  required  by  the  actuator. 
The  dimensionless  specihc  cost  of  transport  is  defined  as: 


,  A  .  l,ol 

cot  = 


Mgd  ’ 


(5.11) 


where  d  is  the  total  distance  travelled  by  the  system,  and  Wact  tot  is  the  total  work  done 
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by  the  actuator  for  the  entire  stance  phase.  To  avoid  power  regeneration,  it  is  appropriate 
to  compute  Wact  tot  as  the  integral  of  the  unsigned  power: 


In  our  case,  the  total  cost  of  transport  can  be  computed  for  one  jump  (apex-to- 
apex),  and  therefore  the  distance  travelled  d  is  equal  to  the  distance  travelled  from 
one  initial  apex  state  to  the  next.  The  work  done  by  the  actuator  is  a  combination  of 
the  work  done  to  cancel  the  nonlinearity  as  in  (5.1),  and  the  work  required  to  move 
the  actuator  at  the  desired  constant  values  £Cii  and  £Ct 2.  Fig.  5.10a  shows  the  cost  of 
transport  computed  for  symmetric  jumps  starting  at  different  initial  conditions,  where 
the  only  actuation  considered  is  the  one  due  to  £ni,  while  £Ci\  and  £c,2  are  set  to  zero. 
For  varying  initial  conditions,  cancelling  the  non-linear  terms  in  equation  implies  a  cost 
of  transport  cot  <  0.21.  One  should  note  that  the  considered  jumps  are,  in  general, 
passively  non-symmetric,  and  therefore  the  cost  of  transport  computed  includes  the  cost 
to  drive  many  non-symmetric  jumps  to  be  symmetric.  Fig.  5.10b  shows  instead  the  cost 
of  transport  for  an  assigned  initial  conditions,  considering  both  the  effect  of  £ni  and  of 
varying  values  of  £c\  and  £&.  Also  in  this  case  the  cost  of  transport  is  cot  <  0.26. 
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Figure  5.10:  Fig.  5.10a  illustrates  the  cost  of  transport  computed  for  symmetric  jumps 
starting  at  different  initial  conditions,  where  the  only  actuation  considered  is  the  one 
due  to  £ni,  while  £c,i  and  lc,2  are  set  to  zero.  Fig.  5.10b  shows  instead  the  cost  of 
transport  for  an  assigned  initial  apex  state,  ya  =  1.3  [m],  xa  =  2  [m/s],  and  9td  =  78 
[deg].  Both  the  effect  of  £ni  and  of  varying  values  of  £c\  and  ic 2  are  considered. 
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5.4  Control  Actions:  Online  Computation  of  Opti¬ 
mal  Parameters 

We  start  by  noting  that,  in  order  to  increase  the  dimension  of  the  reachable  state- 
space  to  a  3-dimensional  space,  we  need  to  add  another  degree  of  freedom:  the  touch-down 
angle  9td,  as  illustrated  in  Fig.  5.4.  In  the  formulation  of  our  strategies,  we  will  then 
consider  three  control  parameters:  the  actuator  values  ic\  and  f?C2,  and  the  touch-down 
angle  9td- 

The  main  advantage  of  using  an  approximation  of  the  stance  phase  versus  its  nu¬ 
merical  solution  is  a  reduction  of  computational  time,  which  increases  the  practicality  of 
performing  online  control  actions.  To  give  an  example,  on  a  representative  pool  of  60,000 
initial  conditions  (apex  state  and  touch-down  angle),  we  computed  the  average  time  to 
simulate  the  stance  phase  using  Matlab’s  function  ode45  versus  an  analytical  approxi¬ 
mation.  The  calculations  were  performed  on  a  Microsoft  Windows  based  computer  (Intel 
Core  i7  eight  core  processor  CPU,  2.80  GHz)  using  Matlab  version  R2012a.  While  the 
average  time  for  ode45  was  0.0259  [s],  the  average  time  for  an  approximate  solution  was 
9.4243  *  10-5  [s]:  a  decrease  in  computation  time  of  over  250  times. 

We  start  from  this  preliminary  remark  to  introduce  our  proposed  control  actions. 

5.4.1  Controlling  the  {y,  x}— state  space 

In  this  subsection  we  will  focus  our  attention  on  controlling  the  height  and  velocity 
at  apex  only,  disregarding  the  forward  position.  Therefore,  the  apex  state  becomes  a 
two-dimensional  vector  s  =  {ya,  xa}.  We  use  a  modified  version  of  the  Matlab  function 
fminsearch  (which  optimizes  constrained  problems  using  the  Nelder-Mead  algorithm). 
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At  any  current  apex  state  sn  =  {yn,  xn},  we  compute  the  values  for  the  touch-down 
angle  9td  and  the  two  actuator  values  ic\  and  £c2  that  minimize  in  one  jump  the  distance 
to  a  desired  apex  sdes  =  {ydes,  xdes}.  The  optimization  problem  is  defined  to  be  con¬ 
strained  due  to  the  bounds  on  the  values  taken  by  the  touch-down  angle  and  the  actuator 
displacement. 

At  each  step,  n,  the  cost  function  to  be  minimized,  J(n),  is  defined  as: 

j f  \  1  nn  /  (jjn+l  ~  Vdes)2  T2(xn+ 1  ~  ^des)2  1  n\ 

J\n)  =  100W - -2 - + - ^ - ,  (5.12) 

which  expresses  the  percentage  distance  from  the  next  apex  state  sn+1  =  {yn+ 1,  xn+i} 
to  the  desired  one,  sdes. 

Note  that  the  average  time  to  solve  the  optimization  problem  has  been  computed  to  be 
0.05  [s],  which  is  much  faster  than  the  average  ballistic  apex-to-apex  time.  In  particular, 
the  shortest  apex-to-apex  time  corresponds  to  a  take-off  and  touch-down  angle  6to  =  90 
[deg]  and  Otd  =  90  [deg].  In  order  to  guarantee  the  apex-to-apex  time  to  be  smaller  than 
0.05  [s],  it  is  required  for  the  apex  height  to  be  yap  >  1.0123t?0- 

5.4.2  Adaptive  control  for  steady-state  locomotion 

Due  to  errors  in  the  approximation  (see  Fig.  5.8a,  5.8b  and  5.8c),  the  touch-down  angle 
and  actuator  values  that  minimize  the  cost  function  in  Equation  (5.12)  may  drive  the 
system  to  an  apex  state  that  differs  from  the  desired  one  by  a  certain  amount.  Therefore, 
the  system  will  converge  to  an  apex  state  that  is  not  the  desired  one.  In  order  to  address 
this  problem,  we  propose  here  a  strategy  to  reduce  such  steady-state  error,  driving  the 
system  closer  to  the  desired  apex  state  over  time.  The  strategy  is  summarized  as  follows. 
Let  us  start  from  an  initial  apex  state,  {r/(0),  i(0)},  and  let  us  assume  we  want  to 

ultimately  reach  the  value  {ydes, o,  ides, o}-  At  each  step  n,  we  define  the  error  between 

63 


Approximation  and  two-element  control 


Chapter  5 


the  actual  and  the  approximate  state  as 


A Un  Vn  Vdesfli 


(5.13) 


(5.14) 


where  y  and  x  are  the  actual  apex  height  and  velocity  of  the  system.  At  the  n  —  th  step, 
we  update  the  desired  value  for  the  next  step  ydes,n+ i  and  XdeS)n+  i  to  be 


Vdes,n+ 1  Vdes,n  &l  Ayn 


The  proportional  gains  o\  and  02  are  chosen  to  be  0  <  <J\  <  1,  0  <  o<i  <  1.  The 
desired  apex  state  is  updated  at  each  step,  until  the  errors  A y  =  0  and  Ax  =  0,  and  the 
system  reaches  an  equilibrium.  Fig.  5.11a,  5.11b,  5.11c  and  5. lid  show,  respectively,  the 
percentage  distance  J(n )  from  the  desired  apex  state  after  n  =  1,  3,  6  and  9  jumps,  on  flat 
terrain,  for  7  =  20  and  cri  =  cr2  =  0.8.  Our  proposed  controller  reduces  the  percentage 
distance  J  after  9  jumps  from  a  maximum  of  about  70%  to  a  maximum  of  about  0.3%. 
Fig.  5.12a,  5.12b,  and  5.12c  show  an  example  of  error  reduction  for  3  different  initial 
conditions:  the  error  converges  to  zero. 
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steps  =6 


0.5  1  1.5  2  2.5  3 

x  [m/s] 


steps  =9 


0.5  1  1.5  2  2.5  3 

x  [m/s] 


(c) 


(d) 


Figure  5.11:  These  plots  show  the  percentage  error  J  after  1  (a),  3  (b),  6  (c)  and  9  (d) 
jumps,  for  7  =  20.  The  x-axis  and  y-axis  represent  the  apex  velocity  and  the  apex 
height,  respectively.  Parameters  a\  and  <72  have  been  chosen  to  be  equal  to  0.8. 
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(a) 


(c) 


Figure  5.12:  Evolution  at  each  jump  of  the  apex  errors  Ay  (5.13)  and  Ax  (5.14),  and 
the  cost  function  J  (5.12)  for  inital  and  desired  apex  state  s  =  {ya,  xa}  =  {2.4,  2.8} 
(a),  {1.8,  1.6}  (b),  and  {1.3,  0.8}  (c).  As  we  can  see,  the  errors  converge  to  zero. 
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5.4.3  Foothold  Placement  Control 

One  of  the  key  problems  in  legged  locomotion  on  different  kinds  of  terrain  is  to 
determine  whether  the  foothold,  i.e. ,  the  point  at  which  the  foot  comes  in  contact  with 
the  terrain,  is  considered  safe.  For  example,  we  can  imagine  the  case  of  a  terrain  where 
only  a  specific  set  of  N  footholds,  i  —  {1, . .  .N},  is  allowed,  and  everything  else 

has  to  be  avoided,  as  studied  in  [15].  We  can  then  plan  a  trajectory  to  follow  based 
on  the  knowledge  we  have  of  the  terrain.  In  particular,  at  each  jump  i  we  want  to  find 
a  sequence  of  control  actions  4n  and  £C2,  to  minimize  the  distance  between  the 
desired  footholds,  x^esi  and  the  actual  landing  of  the  foot,  x? : 

di  =11  x{  -xfdesi  ||2, 

where  the  position  of  the  foot  at  touch-down  is  computed  as: 

xf  =  xa  +  Xa\  ~{ya  +  4  cos  9td )  +  4  sin  0TD. 

V  g 

The  foothold  error  is  dependent  on  the  planning  horizon ,  i.e.,  the  number  of  jumps  we 
can  pre-compute.  For  a  planning  horizon  of  length  N,  we  can  write  the  following  cost 
function  to  minimize: 

N 

4  =  E  4 
2—1 

However,  to  ensure  that  the  hopper  maintains  a  certain  desired  height,  ydes,i,  with  respect 
to  the  terrain,  we  can  modify  the  above  cost  to  be: 

N  N 

Jn  =  'y  ^  di  +  wy  y  ](yt  —  ydes,i)  > 

2—1  2—1 
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where  wy  G  [0,  1]  is  the  associated  weight.  Furthermore,  to  ensure  that  the  optimal 
solution  at  the  TV-th  step  does  not  result  in  a  jump  with  negative  forward  velocity, 
another  element  can  be  added  to  the  cost  function,  and  we  obtain: 

N  N 

Jn  =  ^d2i  +w  -  Vdes,i)2  +  w±(xN  -  %_i)2,  (5.15) 

i=  1  i=l 


with  w±  G  [0,  1]. 

Ideally,  one  would  want  an  infinite  planning  horizon:  N  =  oo.  In  fact,  the  longer 
the  horizon,  the  better  the  performance  of  the  optimized  problem  (e.g.,  see  Fig.  5.13). 
However,  the  horizon  length  affects  the  computation  time  required  to  plan  the  desired 


Figure  5.13:  Example  of  foothold  error  for  100  jumps  for  the  case  with  horizon  N  =  3 
(green  bars)  and  N  =  5  (black  bars).  As  expected,  a  higher  planning  horizon  gives 
better  performances  in  terms  of  error. 


trajectory  online.  This  is  where  having  an  approximation  for  the  stance  phase  becomes 
highly  beneficial:  in  fact,  being  able  to  compute  the  optimal  path  via  approximation 
significantly  reduces  the  computational  time,  and  as  a  consequence,  it  is  possible  to 
extend  the  planning  horizon.  However,  one  should  keep  in  mind  that  the  approximation, 
as  such,  carries  an  error:  there  is  then  a  trade-off  between  horizon  length/computation 
time,  and  foothold  error.  Fig.  5.14a  and  5.14b  show  an  example  of  trajectory  planning  on 
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flat  terrain,  comparing  a  case  in  which  the  optimization  problem  has  been  solved  using 
the  approximate  solution  for  the  stance  phase,  and  a  case  in  which  the  numerical  solution 
has  been  used.  As  expected,  the  minimization  via  numerical  solution  has  a  smaller  error 
but  a  higher  computational  time  compared  to  the  minimization  via  approximation  (on 
the  same  horizon  length  N).  To  be  able  to  perform  an  online  computation  of  the  optimal 
control  parameters,  the  computation  time  needs  to  be  much  smaller  than  the  average 
time  during  flight:  t  0.5  [s].  Then,  the  only  viable  option  for  the  numerical  solution 
case  is  N  =  2,  with  average  foothold  error  over  10,000  jumps  of  3.64  [cm].  For  the 
approximate  solution  case,  instead,  one  can  use  N  =  6,  and  the  mean  foothold  error  is 
only  0.45  [cm],  i.e. ,  8  times  smaller  than  the  numerical  case. 
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Figure  5.14:  These  plots  show  the  mean  foothold  error  (a)  and  the  mean  computation 
time  (b)  over  10,000  jumps,  for  trajectory  planning,  for  different  horizon  lengths.  The 
terrain  has  been  chosen  to  be  flat,  with  footholds  drawn  from  a  uniform  distribution 
on  the  open  interval  (0.34),  4))-  The  blue  dotted  bars  refer  to  the  optimization  using 
the  numerically  computed  solution,  while  the  yellow  bars  refer  to  the  optimization 
computed  using  the  approximate  solution. 
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5.5  Applications 

In  this  section  we  show  the  performance  of  our  controller  to  perturbations  in  terrain. 
In  particular,  we  test  recovery  for  perturbations  on  the  terrain  height,  and  we  show  an 
example  of  running  on  rough  terrain. 

The  parameters  and  initial  conditions  used  in  our  simulations  are  defined  in  Table  6.1, 
and  are  based  on  biological  data  for  a  typical  human.  In  [3],  the  relative  spring  stiffness 
7  was  found  to  be  very  similar  for  different  gaits,  such  as  running  and  hopping,  of  various 
animals.  In  particular,  values  were  computed  to  be  between  7.1  and  14.6  for  runners, 
and  7.7  and  13.6  for  hoppers.  Because  of  the  limitations  shown  in  Fig.  5.6a  and  5.6b, 
for  our  simulations  we  use  7  =  20,  which  is  slightly  higher  than  the  average  biological 
value,  but  it  also  matches  a  hardware  prototype  currently  under  development  in  our  lab. 
During  simulations,  if  either  the  total  actuator  displacement  £act  =  £ni  +  £c  or  the  total 
actuator  velocity  vact  =  vni+vc  required  were  exceeding  the  maximum  values  allowed,  the 
actuator  was  assumed  to  saturate  its  maximum  allowed  value  for  £act  or  vact ,  respectively. 
Furthermore,  to  acknowledge  the  time  to  solve  the  optimization  problem,  we  limit  our 
simulations  to  initial  apex  heights  and  touch  down  angles  that  corresponds  to  a  time 
during  flight  tf  >  0.15  seconds. 

5.5.1  Recovery  from  perturbations 

We  consider  the  set  of  initial  conditions  in  Table  6.1.  We  test  the  recovery  capabilities 
of  our  controller  when  the  active  SLIP  encounters  an  unexpected  (positive  or  negative) 
perturbation  on  the  terrain  height  of  up  to  50%  of  the  leg  length  £0. 

At  each  apex  state,  we  use  the  strategy  in  5.4.1  to  compute  the  optimal  values  for 
£ci,  iC2  and  dTD  for  flat  terrain.  Once  the  leg  touches  the  ground  with  the  computed 
touch-down  angle  and  the  desired  £&,  the  strategy  in  5.4.1  is  simulated  again  during 
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Simulation  Parameters 


9  = 

9.81  m/s2 

M  = 

80  kg 

4  = 

1  m 

4,0  = 

0.5  m 

^ k,min 

0.05  m 

^ act  ^ 

[-0.1,  0.1]  m 

4  6 

[-0.05,  0.05]  m 

Vact  < 

1  m/s 

Vc  = 

0.5  m/s 

y  e 

[1,  2.5]  m 

x  G 

[0.5,  3]  m/s 

7  = 

20 

Table  5.1: 


the  first  half  of  the  stance  phase,  this  time  to  compute  only  ic 2  to  take  into  account 
the  encountered  perturbation  on  the  terrain  height.  Note  that  the  search  of  this  second 
value  takes  on  average  0.01  [s],  which  is  much  smaller  than  the  average  time  required 
for  the  first  half  of  the  stance  phase,  and  therefore  it  can  be  realistically  implemented. 
Fig.  5.15a,  5.15b,  5.15c,  and  5.15d  show  the  number  of  jumps  necessary  for  the  system 
in  order  to  return  within  1%  of  the  initial  apex  state,  in  the  case  of  positive  or  negative 
perturbations.  Our  controller  is  robust  to  perturbations  of  varying  magnitudes,  with 
ability  to  recover  in  up  to  8  jumps.  Clearly,  these  results  are  not  only  due  to  the  control 
strategy,  but  are  partly  affected  by  actuator  limits.  One  can  expect  the  number  of 
jumps  to  recover  to  increase  or  decrease  if  the  actuator  limits  are  more  or  less  stringent. 

Fig.  5.16a  and  5.16b  show  an  example  of  apex  recovery  for  an  unforeseen  drop  of 
magnitude  50%  of  the  leg  length. 
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step  height  =-0.2  step  height  =-0.5 


x  [m/s]  x  [m/s] 

(c)  (d) 


Figure  5.15:  Number  of  jumps  (colorbar)  to  reach  1%  of  desired  value,  for  terrain 
perturbation  of  magnitude  (a)  —0.2  [nr],  (b)  —0.5  [nr],  (c)  0.2  [nr],  (d)  0.5  [nr].  Simu¬ 
lation  parameters  are  chosen  as  per  Table  6.1.  Note  that  in  case  of  a  positive  terrain 
perturbation  (subplots  (c)  and  (d)),  the  apex  height  y  has  been  chosen  to  leave  enough 
room  to  the  leg  to  swing  during  flight  without  colliding  with  the  terrain. 
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step  # 
(a) 


Figure  5.16:  Plot  (a)  shows  the  trajectory  of  the  active  SLIP  for  initial  apex  state 
y  =  2  [m],  x  =  1.3  [m/s],  and  negative  terrain  perturbation  of  magnitude  0.6  [m].  The 
blue  dotted  line  represent  the  trajectory  of  the  mass,  the  purple  circles  the  desired 
apex  height.  The  green  solid  line  is  the  terrain  height  and  the  red  x  marks  are  the 
landing  points  of  the  foot.  Plot  (b)  shows  the  percentage  error  of  apex  height  (red 
circle)  and  velocity  (blue  diamond)  at  each  jump  after  the  drop.  We  can  see  that  the 
controller  reaches  and  remains  within  1%  of  the  original  apex  state  in  5  jumps. 

5.5.2  Hopping  on  rough  terrain 

In  this  section  we  show  an  example  of  the  active  SLIP  model  hopping  on  rough  terrain. 

Additionally,  we  assume  that  the  estimates  of  upcoming  terrain  height  are  faulty,  and  we 
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want  to  maintain  the  same  forward  velocity  and  the  same  distance  from  the  terrain  with 
respect  to  the  last  jump.  As  we  can  see  from  Fig.  5.17a,  5.17b,  and  5.17c,  the  system  is 
able  to  successfully  hop  on  a  rough  terrain  with  a  maximum  magnitude  landing  height  of 
0.26  [m]  and  a  maximum  perturbation  of  0.46  [m],  i.e. ,  46%  of  the  leg  length.  Consistent 
with  what  is  in  Fig.  5.8b  and  5.8c,  the  error  on  x  is  on  average  higher  than  the  error  on 
y,  and  in  this  example  they  both  do  not  exceed  4%.  Note  that,  since  the  terrain  height 
varies  continuously,  the  energy  of  the  system  varies  at  each  jump. 

The  strategy  of  recomputing  £c2  during  the  first  portion  of  the  stance  phase  after 
encountering  a  perturbation  can  dramatically  improve  the  performance  of  our  controller. 
Indeed,  performance  on  the  same  terrain  has  been  studied  using  the  computation  during 
flight  only  versus  with  controller  update  via  recomputing  of  £&  during  the  first  half  of 
stance.  Fig.  5.17c  shows  that  the  cost  function  J  (computed  as  in  (5.12))  is  significantly 
smaller  than  the  case  without  controller  update. 
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Figure  5.17:  Plot  (a)  shows  the  trajectory  of  the  active  SLIP  hopping  on  a  random¬ 
generated  rough  terrain.  The  controller  acts  to  maintain  a  constant  apex  height  with 
respect  to  the  terrain  of  y  =  1.5  [m],  and  a  constant  forward  velocity  of  x  =  2  [m/s]. 
The  blue  dotted  line  represent  the  trajectory  of  the  mass,  the  purple  circles  the  desired 
apex  height.  The  black  dotted  line  is  the  expected  terrain  height,  while  the  green  solid 
line  is  the  actual  terrain  height  and  the  red  x  marks  are  the  landing  points  of  the  foot. 
Plot  (b)  shows  the  percentage  error  of  apex  height  PEy  (circle)  and  velocity  PE± 
(diamond)  at  each  jump.  Plot  (c)  shows  the  cost  J  as  in  (5.12)  for  the  case  without 
(blue  star)  and  with  (red  square)  controller  update  during  the  first  part  of  the  stance 
phase. 
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5.5.3  Landing  on  feasible  footholds 

We  show  here  an  example  of  hopping  on  a  set  of  predefined  footholds. 

The  set  of  footholds  on  an  uneven  terrain  have  been  generated  taking  points  from  a 
uniform  distribution  between  0.6  [m]  and  1.2  [m].  The  cost  function  the  be  minimized 
is  defined  in  Eq.  (5.15),  with  weights  wy  =  0.8  and  w±  =  1.  We  used  a  least-square 
algorithm  to  find  the  optimal  solution  for  a  planning  horizon  of  N  =  4  steps,  for  a  total 
of  50  steps.  Note  that,  according  to  the  specifics  of  our  minimization  algorithm,  the 
hopper  was  allowed  to  skip  some  footholds  if  pertinent.  Fig.  5.18a  shows  a  close-up  of 
the  output  trajectory,  while  Fig.  5.18b  shows  at  each  step  the  foothold  error  (i.e.,  the 
distance  from  the  desired  to  the  actual  foothold),  the  forward  velocity  and  the  apex 
height.  The  desired  apex  height  with  respect  to  the  terrain  was  set  as  ydes  —  1-4  [m].  As 
we  can  see,  over  50  steps,  the  mean  foothold  error  was  2.05  [cm]. 
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(a) 


Mean  foothold  error  =  2.05  [cm] 


20  25  30  35  40 


step  # 


(b) 

Figure  5.18:  Plot  (a)  shows  the  trajectory  of  the  point  mass  (blue  dotted  line)  on 
a  terrain  (green  line).  The  white  circles  are  the  desired  footholds,  while  the  red  x 
marks  are  the  actual  placements  of  the  foot.  As  noted  in  the  text,  our  minimization 
algorithm  allows  for  some  footholds  to  be  skipped.  Plot  (b)  from  top  to  bottom  shows 
the  foothold  error  (red  bars),  the  forward  velocity  (yellow  bars)  and  the  apex  height 
with  respect  to  the  terrain  height  (green  bars).  The  dotted  blue  line  represents  the 
desired  apex  height. 
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5.5.4  Extension  to  general  planar  hopper  models 

The  actuated  SLIP  model,  while  encompassing  the  general  behavior  of  hopping  spring- 
mass  systems,  is  a  highly  idealized  model.  It  is  then  logical  to  ask  how  the  strategy 
presented  in  this  work  would  fare  when  applied  to  models  that  match  more  accurately 
real  hardware  prototypes.  To  answer  this  question,  the  two-step  strategy  has  been  tested 
on  a  more  realistic  hopper  model,  shown  in  Fig.  5.19,  similar  to  the  hopper  studied  in  [10]. 
The  leg  has  mass  mi  and  moment  of  inertia,  .//,  while  the  body  has  mass  and  inertia  M 
and  J.  Body  and  leg  are  connected  by  a  joint  at  the  hip,  where  an  actuator  can  apply  a 
torque,  r/wp.  Equations  of  motion  for  this  system  can  be  easily  computed,  for  example  via 


Figure  5.19:  Planar  hopper  model.  The  body  has  mass  M  and  moment  of  inertia  J, 


and  center  of  mass  at  the  hip.  <j)  is  the  body  angle.  The  leg  has  unsprung  mass  mi  and 
moment  of  inertia  J;,  and  its  center  of  mass  is  located  at  distance  ti  from  the  foot.  £ 
is  the  total  leg  length,  and  9  the  angle  that  the  leg  forms  with  respect  to  the  ground 
during  contact.  The  spring  has  stiffness  k.  lact  denotes  the  length  of  the  series  elastic 
actuator.  An  actuator  at  the  hip  can  apply  a  torque  Thip  between  the  body  and  the 
leg. 

the  Lagrangian  method.  During  flight,  a  PD  controller  is  applied  to  the  the  hip  motor 
to  position  the  leg  at  the  desired  angle,  and  the  center  of  mass  of  the  system  follows  a 
ballistic  trajectory.  The  effect  of  the  energy  loss  at  impact  on  the  velocities  of  the  leg 
length  and  leg  angle  is  approximated  to  determine  the  initial  (post-impact)  conditions 
of  the  stance  phase.  During  stance,  the  rotation  of  the  body  affects  the  dynamics  of  the 
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leg  angle,  while  the  dynamics  of  the  leg  length  is  not  affected  and  can  be  solved  as  in 
Section  5.1.  To  keep  a  body  attitude,  i.e.,  to  prevent  the  body  from  tipping  forward  or 
backward  during  stance,  a  torque  is  applied  at  the  hip.  Therefore,  an  approximation 
of  the  equation  of  motion  for  the  angular  displacement  over  time,  6{t),  is  given  by  the 
approximation  computed  in  Section  5.1,  with  the  addition  of  a  term  that  approximates 
the  effect  of  the  motion  of  </>(£).  Note  that  this  is  still  an  approximation,  and  we  exploit 
the  actuator  at  the  hip  by  applying  an  additional  torque  to  control  the  system  to  follow 
the  approximated  trajectory.  The  series  elastic  actuator  is  modelled  as  a  function  of  input 
current,  with  a  damping  and  friction  term.  Details  of  the  approximation,  controllers  and 
parameter  values  are  illustrated  in  Appendix  A. 

The  control  strategy  has  been  tested  on  a  flat  terrain.  Every  10  jumps  the  target 
apex  height  and  velocity  have  been  changed  to  random  values  taken  form  a  uniform 
distribution  between  ±0.1  [m]  and  ±0.3  [m/s]  from  the  previous  target  state,  respectively. 
Fig.  5.20a  and  5.20b  show  the  results  for  100  jumps.  The  system  is  controlled  to  reach 
simultaneously  the  apex  height  and  forward  velocity  desired. 
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Figure  5.20:  Data  showing  the  apex  state  reached  for  a  set  of  100  jumps.  Plot  (a)  and 
(b)  show,  respectively,  the  apex  height  and  forward  velocity  reached  (red  dots)  and 
the  desired  values  (black  line). 
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5.6  Conclusions 

In  this  chapter  we  proposed  a  strategy  for  actuator  displacement  of  the  active  SLIP 
model.  On  one  hand,  it  allows  us  to  analytically  solve  the  equation  that  describes  the 
leg-length  dynamics  during  stance  via  partial  feedback  linearization.  On  the  other,  by 
dividing  the  stance  phase  in  two  parts  defined  by  the  point  of  maximum  leg  compression 
and  setting  two  different  actuation  values  for  each  part,  it  allows  us  to  add  or  remove 
energy  and  to  modify  the  upcoming  apex  state  to  span  an  open  set  within  the  reachable 
apex  state  set.  Our  strategy  was  tested  to  deal  with  terrain  perturbations,  and,  for  a 
set  of  system  parameters,  we  quantitatively  define  the  number  of  jumps  necessary  for 
full  recovery.  We  additionally  extend  the  proposed  approximation  to  the  case  of  a  more 
realistic  hopper  model  with  leg  mass  and  body  inertia. 
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As  seen  in  the  previous  chapters,  actuating  the  leg  has  the  effect  of  modifying  the  apex 
state  that  would  otherwise  be  reached  passively.  In  particular,  each  specific  actuation 
policy  differently  modifies  the  reachable  space  of  the  active  SLIP  model.  In  this  chapter, 
we  want  to  answer  the  following  questions.  “Is  it  possible  to  characterize  the  effect  that 
the  actuation  has  on  the  following  apex  state?  How  does  applying  different  amounts  of 
actuation  at  various  instances  during  the  stance  phase  affect  the  apex  state?  Furthermore, 
is  there  an  actuation  strategy  that  maximizes  the  reachability  space?” 

Displacing  the  actuator  during  the  stance  phase  affects  all  three  dimensions  of  the 
next  apex  state  reached.  However,  any  leg  actuation  strategy  applied  to  its  motion 
can  only  control  two  of  the  three  dimensions  simultaneously,  as  shown  in  Chapter  3. 
For  this  reason,  here  we  choose  to  consider  the  apex  height  and  forward  velocity  only, 
disregarding  the  forward  position  x,  which  is  useful  mostly  in  path  planning  scenarios. 
Hence,  we  will  define  the  reduced  apex  state  as  S  =  {ya,  xa}.  In  this  work,  we  are 
interested  in  determining  how  the  actuator’s  displacement  affects  the  reachable  space. 
Our  study  does  not  aim  to  propose  a  leg-placement  strategy,  but  rather  to  understand 
the  effect  that  the  actuation  has  during  the  stance  phase.  Therefore,  we  limit  our  analysis 
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to  a  reduced  reachable  space,  restricted  to  the  case  where  the  touch-down  angle,  0TD,  is 
fixed.  The  reachable  space  defined  in  (3.4)  reduces  to: 


T^-{Sq,  0td)  {Si  |  3  £act{t)  ■  X{Sq,  Otd-,  £a{t))  —  *S'i}- 


Obviously,  H(So,9td)  Q  7£(so)- 

As  mentioned  in  Chapter  3,  adding  the  series  actuator  has  the  effect  of  extending  the 
reachable  set  to  a  2-d  surface,  whose  shape  is  determined  by  the  shape  of  the  function 
(■act  (t)  during  the  stance  phase.  A  few  leg- actuation  functions  have  been  proposed  and 
can  be  found  in  the  literature.  Among  them,  we  choose  two  functions  that  we  consider 
relevant  because  they  provide  an  insightful  representation  of  how  different  leg-actuation 
motions  can  affect  the  locus  of  all  the  points  that  can  be  reached  in  one  step,  indeed, 
we  want  here  to  compare  their  reachability,  to  determine  the  best  course  of  action  when 
choosing  an  actuation  strategy. 

In  [13],  J.  Schmitt  and  J.  Clark  propose  a  leg-length  actuation  inspired  by  studies  on 
the  net  energy  production  of  muscles  groups  in  humans: 

£act(t)  Cct, 0  £ dev  Sm  (cut),  (6.1) 

where  £act,o  is  the  initial  actuator  length,  £dev  is  the  amplitude  of  the  variation  of  the 
actuator  length  from  its  nominal  length,  and  u  is  the  frequency  of  the  actuator’s  motion. 
t  is  the  time  from  the  beginning  of  the  stance  phase,  and  t  —  0  corresponds  to  the  touch¬ 
down  time.  Let  us  use  7 Zs(S0,  Otd )  to  indicate  the  reachable  space  computed  with  the 
actuator  movement  described  in  (6.1). 

In  Chapter  5,  the  stance  phase  is  divided  in  two  parts  marked  by  the  point  of  maximal 
spring  compression  (which  is  equivalent  to  £{t)  =  0).  The  actuator  is  then  moved  from 
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its  initial  position  to  a  desired  value  4,i  during  the  first  half  of  the  stance  phase,  and  a 
desired  value  4,2  during  the  second  half: 


Lctif) 


j  4,i,  Vt  :  £(t)  <  0 
^4,2,  Vt  :  £(t)  >  0. 


(6.2) 


According  to  the  finite  response  of  a  physical  actuator,  the  actuator  does  not  move 
instantaneously,  but  it  is  assumed  to  move  with  its  maximum  allowable  velocity  and 
acceleration.  Let  7Zp(S0,  0TD)  be  the  reachable  space  computed  with  the  actuator  move¬ 
ment  (6.2),  as  described  in  Chapter  5. 

How  do  the  reachable  spaces  7Zs(S0j9td)  and  7Zp(S0,6td )  compare?  And,  how  do 
they  compare  with  respect  to  all  the  possible  apex  states  reachable  in  one  step  for  other 
actuator’s  motions?  To  answer  to  these  questions,  we  compare  the  two  reachable  sets 
with  the  set  of  all  reachable  states  obtained  by  moving  the  actuator  at  any  time  during 
the  stance  phase  in  any  direction.  Let  us  call  this  space  TZtot(So,  &td)-  Fig-  6.1  shows 
the  reachability  space  of  the  three  methods  considered,  starting  from  the  same  initial 
conditions  at  touch-down,  and  with  the  same  actuator’s  limitation.  The  actuator  has 
been  modelled  as  in  Equation  (3.5).  As  we  can  see,  the  area  of  7 Ztot  is  bigger  than  the 
areas  of  7 Zs  and  7 Zp.  Furthermore,  the  convex  ellipsoid-like  shape  of  7 Ztot  presents  the 
advantage  of  reaching  any  point  close  to  a  desired  state.  This  suggests  that,  in  order  to 
reach  a  wider  range  of  apex  states,  it  is  beneficial  to  consider  a  strategy  that  involves 
updating  the  actuator  movement  throughout  the  entire  stance  phase,  instead  of  choosing 
a  priori  some  fixed  parameters.  Additionally,  it  was  shown  in  Fig.  3.2  that  the  complete 
reachable  space  is  a  function  of  the  actuator’s  characteristics.  As  a  result,  in  the  strategy 
proposed  in  this  work,  we  will  focus  on  moving  the  actuator  with  its  maximum  allowed 
parameters. 
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[fotatiSo,  @Td) 
\1Zp(Sq,  Otd) 
}Ks(Sq,  Otd) 


4  4.2 

x[m/s\ 


Figure  6.1:  Reachable  space  of  a  passively  symmetric  jump  7Z(Sq,  Otd )  computed  for 
So  =  (1.3  [m],  4  [m/s]),  Otd  =  122.14  [deg],  leg  length  £q  =  1  [m],  and  7  =  10.  Max¬ 
imum  actuator  velocity  allowed  vmax  =  .5  [m/s]  and  maximum  acceleration  constant 
kacc  =  10  [m/s2],  respectively.  Defining  Atot ,  Ap  and  As  the  areas  of,  respectively, 
tot,  Tip,  and  TZS,  we  have  that  Ap  =  0.61{L4tot,  and  As  =  0.445^4tot.  Note  that,  while 
p  C  ittot,  this  does  not  necessarily  hold  true  for  7ZS.  since  its  actuator  dynamics  (6.1) 
start  with  nonzero  velocity. 
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6.1  Actuator  movement 


In  this  section,  we  aim  to  understand  and  characterize  the  effects  that  actuation 
motion  has  on  the  next  apex  state  reached.  We  will  show  the  relationship  between 
moving  the  actuator  at  different  times  during  the  stance  phase  and  the  variation  of  the 
apex  state  thereby  reached. 

Our  objective  is  to  find  an  actuator  movement  strategy  to  (i)  move  the  system  to 
a  desired  apex  and  (ii)  counteract  the  effects  of  perturbations  (on  terrain  height  or  on 
initial  position  at  apex).  As  previously  stated  in  Section  2.2,  compressing  or  extending  the 
actuator  during  the  stance  phase  corresponds  to  an  increase  or  decrease  of  the  system’s 
energy:  the  additional  potential  energy  stored  in  the  spring  through  actuation  will  then 
be  transformed  into  additional  kinetic  energy  throughout  the  stance  phase,  until  take-off, 
and  vice  versa.  However,  the  system’s  energy  variation  is  a  function  of  the  values  of  the 
states  of  the  system  when  actuation  is  performed.  Therefore,  to  drive  the  system  from 
an  initial  state  S0  =  { y0 ,  i0}  to  a  desired  apex  state,  Sdes  =  {Vdes,  Xdes},  it  is  crucial  to 
understand  the  effects  that  the  particular  time  during  stance  when  actuation  happens 
has  on  the  evolution  of  the  states,  and,  therefore,  on  the  position  of  the  reached  apex. 

The  lack  of  closed-form  solution  for  the  system’s  dynamics  constrains  us  to  perform 
a  numerical  study.  Toward  generality  in  this  particular  analysis,  we  consider  a  step-type 
actuation,  i.e. ,  we  assume  that  the  actuator  can  be  instantaneously  moved  to  reach  a 
specified  value,  £act,des- 


4  ct(t) 


0,  Vt  <  ts 

£ act,des’)  ts 


(6.3) 


where  ts  is  the  time  at  which  the  step  actuation  occurs.  From  an  inital  apex  state 
and  touch-down  angle,  we  computed  the  passively-reached  state,  i.e.,  the  state  reached 
without  actuation  Spass  =  { ypaSs ,  xpass } .  Given  a  positive  and  negative  step  actuation, 
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we  computed  the  apex  reached  as  a  function  of  the  time  at  which  the  positive  or  negative 
actuation  was  applied:  Sr(£act}des,ts)  =  {y(Lct,des,  ,ts),  i(£acttdes,  ts)}.  Fig.  6.2  shows  an 
example  of  the  variation  of  the  two  components  of  the  reached  apex  with  respect  to  the 
passively-reached  apex,  as  a  function  of  the  time  ts  at  which  the  step  actuation  is  applied: 


£-^yr{,£act,desi^s)  (jjr  {£-act,des  ?  A  )  Vpass  )/4, 

(  £-act,des  i  ts)  ^r{^£-act,desi^s)  %  pass')  /  ‘\j , fj  ^  0  • 

where  the  states  are  normalized  to  be  dimensionless  for  comparison  purposes.  As  we  can 
see,  the  variation  of  yr  follows  the  sign  of  the  step  input:  any  positive  actuation  will  result 
in  an  increase  of  the  next  apex  height,  whereas  any  negative  actuation  will  result  in  its 
decrease,  regardless  of  the  time  at  which  the  actuation  is  performed.  Dissimilarly,  the 
evolution  of  xr  initially  increases  or  decreases  according  to  the  sign  of  the  actuation,  and 


Figure  6.2:  This  figure  shows  the  variation  of  A yr  and  Axr  as  a  function  of  the  time 
ts  at  which  the  step  actuation  is  applied  to  the  system.  The  blue  lines  refer  to  a 
positive  actuation  £act,des  =  0-3  [m],  while  the  red  lines  refer  to  a  negative  actuation 
£act,des  =  —0.3 [to],  for  a  model  with  Iq  =  1  [m]  and  7  =  20.  The  vertical  dotted  line 
signals  the  time  at  which  the  spring  reaches  its  maximum  compression. 
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then  increases  of  decreases  opposite  of  the  sign  of  the  actuation.  By  lacking  analytical 
solution  for  the  stance  phase,  it  has  not  been  possible  to  pinpoint  the  exact  time  at  which 
the  variation  of  apex  state  happens.  However,  we  can  see  that  the  change  in  direction 
for  Axr  happens  when  the  actuator  moves  close  to  the  time  corresponding  to  the  point 
of  half  stance  with  respect  to  the  spring  movement,  i.e. ,  when  the  spring  reaches  its 
maximum  compression.  We  can  then  summarize  the  above  observations  as: 


•  First  half: 

o 

A 

o 

53 

* 

y  T,  and 

*  tact  <  0: 

y  i,  and 

xt, 

•  Second  half: 

*  tact  >  0: 

V  T,  and 

X  t, 

*  tact  <  0: 

y  J,,  and 

X  J,. 

Another  important  insight  we  can  gather  from  Fig.  6.2  concerns  the  different  rates  at 
which  y  and  x  vary  throughout  stance.  While  providing  actuation  during  the  first  half  of 
stance  has  the  biggest  effect  on  the  variation  of  the  x  component  of  the  next  reached  apex 
state,  the  second  half  has  little  effect  on  it.  On  the  contrary,  the  y  component  is  influenced 
throughout  the  entire  stance  phase,  with  maximum  variation  during  the  second  part. 
This  is  an  indicator  and  a  confirmation  of  our  idea  that  varying  the  actuator  throughout 
stance  gives  the  greatest  control  authority  to  reach  a  desired  state,  compared  to  choosing 
any  fixed  time  during  stance  at  which  to  activate  the  actuation.  This,  for  example, 
also  explains  and  is  confirmed  by  the  bowtie-like  shape  of  the  reachability  set  1ZS  in 
Fig.  6.1:  in  order  to  modify  y  without  much  variation  of  the  otherwise  passively  reached 
value  of  x,  the  actuation  should  mostly  be  applied  during  the  first  half  of  the  stance 
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phase.  However,  the  choice  of  the  sinusoidal  function  (6.1)  with  constant  amplitude  and 
frequency,  distributes  the  actuation  almost  evenly  throughout  the  entire  stance  phase, 
precluding  the  possibility  of  favoring  one  part  of  stance  with  respect  to  another. 

These  observations  consider  only  one-step  actuator  movement  during  stance,  but 
can  be  generalized  to  the  case  of  subsequent  actuation  movements,  with  more  realistic 
dynamics,  as  for  example  (3.5).  We  accomplish  this  by  updating  the  actuator  value 
every  few  time-steps  during  the  stance  phase.  Let  us  define  t  as  the  time  elapsed  from 
the  beginning  of  the  stance  phase  (where  to  =  0  is  the  touch-down  time),  and  divide 
the  stance  phase  in  time  intervals  of  length  St.  Let  Sn  be  the  apex  state  that  would  be 
reached  if  the  actuator  were  to  be  moved  as  a  certain  function  fn(t)  from  touch-down  to 
time  nSt,  and  then  stopped  at  the  current  value  for  the  rest  of  the  stance  phase: 


Sn 


-act 


j  fn(t), 

[  fn{nSt), 


if  t  &  [to,  nSt] 
if  t  >  nSt 


Then,  Sn+ 1  will  be  the  apex  state  reached  if  the  actuator  were  to  be  moved  as: 


I  fn+i(t),  if  t  e  [t0,  (n  +  l)St] 

^n+l  ^  tact  \ 

|yn+i((n+  l)St),  if  t  >  (n  +  l)St 
where  fn+i(t)  =  fn(t),  for  t  G  [0,  nSt}. 

Fig.  6.3a  and  6.3c  show  the  position  of  the  next  apex,  Sn,  when  moving  the  actuator 
as  respectively  in  Fig.  6.3b  and  6.3d  for  values  of  n  =  1,2,...  until  take-off.  A  time 
interval  of  length  St  =  0.02  [s]  has  been  chosen.  As  we  can  see,  and  according  to  earlier 
observations,  at  each  time  step  n,  compressing  the  spring  results  in  an  increase  of  yn, 
while  extending  the  spring  results  in  a  decrease  of  yn  with  respect  to  its  value  at  the 
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(n  —  1)— th  step.  Conversely,  any  spring  compression  translates  in  a  decrease  or  increase 
of  xn  if  this  happens,  respectively,  during  the  first  or  second  half  of  the  stance  phase.  Vice 
versa,  an  extension  of  the  spring  applied  during  the  first  or  second  half  of  the  stance  phase 
results  in  an  increase  or  decrease,  respectively,  of  xn.  This  is  summarized  graphically  in 
Fig.  6.4. 
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(a)  (b) 


Figure  6.3:  Subfig.  6.3a  and  Subfig.  6.3c  show  the  next  apex  states  Sn,  n  =  1,2,... 
until  take-off,  for  the  different  actuation  motions  shown  in  Subfig.  6.3b  and  6.3d.  The 
blue  line  and  the  dotted  red  line  in  Subfig.  6.3a  and  6.3c  are  the  apex  states  reachable 
if  the  actuator  motion  follows  the  solid  blue  and  dotted  red  trajectories  in  Subfig.  6.3b 
and  6.3d,  respectively.  The  black  star  symbol  represents  So,  i.e.,  the  passively  reached 
state.  The  yellow  diamonds  represent  the  apex  states  reached  following  the  dotted 
red  actuator  trajectory  in  Subfig.  6.3b  and  6.3d,  while  the  green  squares  represent  the 
apex  states  reached  if  the  actuator  follows  the  solid  blue  trajectory  in  Subfig.  6.3b 
and  6.3d. 
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4  First  half  stance 


y 


X 


4  Second  half  stance 


X 


i — i Lct(ri5t)  >  £act((n  ~  1  )St) 

'---'£act(n8t)  <  £act{(.n  -  1  )8t) 

Figure  6.4:  Position  of  the  apex  state  that  would  be  reached  providing  a  positive  or 
negative  actuation  during  the  following  interval  [(n  —  l)5t,  n5t ].  Note  that  the  slope  of 
the  directional  curve  from  Sn- ±  to  Sn  is  a  function  of  the  particular  actuation  motion, 
the  system’s  parameters  and  the  initial  condition  Sn- ±. 
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6.2  Control  Strategy 

We  introduce  here  our  control  strategy  to  reduce  the  error  between  the  reached  and 
a  desired  apex  state,  first  giving  an  example,  and  then  explaining  in  detail  our  proposed 
algorithm. 

The  objective  of  our  control  strategy  can  be  summarized  as  follows. 

Let  So  =  { 2/0 ,  io}  be  the  initial  apex  state,  and  Otd  the  angle  of  the  leg  at  touch¬ 
down.  Let  us  define  Sdes  —  {lldes^des}  as  the  desired  apex  state,  and  Sr  =  {yr,  xr }  the 
apex  state  reached  with  a  certain  actuation  motion.  Then,  we  want  to  find  an  actuation 
motion  strategy  to  reduce  the  Euclidean  distance  between  the  desired  apex  and  the 
actual  apex  reached.  In  order  to  compute  the  Euclidean  distance  over  variables  with 
different  measurement  units  (apex  height  and  velocity),  we  reduce  the  states  to  their 
dimensionless  counterpart.  By  defining  the  dimensionless  time  r  =  \Jg/£ot,  we  can  write 
the  dimensionless  apex  height  as  y  =  y /i 0,  and  the  dimensionless  apex  forward  velocity 
as  x  =  x/yfgC o-  The  dimensionless  apex  state  is  S  =  {y,  x}.  The  error  is  then  computed 
as: 

Err  =||  Sdes  -  Sr  ||2  •  (6.4) 

We  explain  here  our  proposed  strategy  to  minimize  error  during  stance. 

During  the  stance  phase,  let  us  define  each  time  interval  of  size  5t  as  A tn  =  [n5t,  (n  + 
l)5i],  n  —  1,2,...  until  take-off.  Our  control  strategy  for  actuator  movement  is  to  update 
the  actuator’s  value  every  time-step  5t  of  the  stance  phase,  according  to  predictions  of 
next  apex  state  based  on  the  current  state.  This  means  that  at  each  time  interval 
A tn,  a  new  value  for  the  actuator  displacement  at  the  next  time  interval,  £act(Atn+i),  is 
computed. 

At  each  interval  A tn,  we  can  compute  the  apex  state,  Sn,  that  the  SLIP  model  would 
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reach  if  the  actuator  remains  at  the  current  value  (computed  at  A tn_i): 

Lct{t)  =  £act((n-l)St),  Vt  >  n5t, 

i.e.,  the  apex  state  reachable  without  further  actuator  movement.  This  is  motivated  by 
the  observations  highlighted  in  Section  6.1,  i.e.,  that  activating  the  actuator  during  the 
first  part  of  the  stance  phase  has  the  effect  of  increasing  or  decreasing  the  value  of  the 
forward  height  at  the  next  apex,  y,  while  at  the  same  time  decreasing  or  increasing, 
respectively,  the  value  of  the  next  apex’s  forward  velocity,  x.  By  contrast,  moving  the 
actuator  during  the  second  half  of  the  stance  phase  has  the  effect  of  increasing  or  de¬ 
creasing  the  values  of  both  y  and  x  at  the  next  apex  state.  Hence,  to  determine  how  to 
control  the  actuator’s  displacement  at  the  next  time  step,  it  is  imperative  to  determine 
the  position  of  Sn  with  respect  to  the  position  of  the  desired  apex  state  Sdes •  At  the  same 
time,  the  choice  of  the  next  value  for  the  actuator,  £act(Atn+\) ,  is  computed  predicting 
the  position  of  Sn+ ±  when  moving  the  actuator  for  one  step  At  to  one  direction  or  the 
other,  i.e.,  compressing  or  extending  the  spring. 

For  example,  let  us  assume  that  the  relative  position  of  Sdes  with  respect  to  Spass  is 
given  in  Fig.  6.5  /.  From  what  is  graphically  summarized  in  Fig.  6.2  and  Fig.  6.4,  we  can 
already  expect  that  the  actuator  motion  necessary  to  reach  the  desired  state  will  consist 
of  either  an  initial  spring  extension  or  compression,  followed  by  a  compression  during 
the  second  half  of  the  stance  phase.  What  we  aim  to  accomplish  with  our  algorithm, 
is  to  control  the  amount  of  the  desired  compression/extension,  by  updating  its  value 
throughout  the  entire  stance  phase.  At  each  time-step  during  the  first  half  of  the  stance 
phase,  we  will  compute  the  next  apex  reached  if  the  actuator  were  to  be  compressed  or 
extended  during  the  next  time-step,  Fig.  6.5  II,  kept  at  the  reached  value  for  the  rest  of 
the  first  half  of  the  stance  phase,  and  subsequently  compressed  during  the  entire  second 
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half  of  the  stance  phase,  as  in  Fig.  6.5  III .  By  computing  and  averaging  the  distance 
of  the  resulting  apex  states  to  the  desired  state,  we  can  compute  how  much  to  move  the 
actuator  at  the  next  step,  Fig.  6.5  IV.  Once  the  dynamics  enter  into  the  second  half 
of  the  stance  phase,  Fig.  6.6  /,  we  can  start  compressing  the  actuator:  we  compute  the 
states  the  system  would  reach  if  the  actuator  were  to  move  at  its  maximum  and  at  half 
its  capabilities  (e.g.,  actuator  acceleration)  at  the  next  time  step,  and  then  kept  at  the 
current  value  for  the  rest  of  the  stance  phase,  Fig.  6.6  II .  Based  on  an  average  of  the 
distance  of  these  two  predicted  states  with  respect  to  the  desired  apex,  the  actuator’s 
action  to  implement  at  the  following  time  step  can  be  computed,  Fig.  6.6  III.  Final 
result  of  our  algorithm  is  shown  in  Fig.  6.6  IV. 

Since  at  each  step  the  algorithm  predicts  where  the  next  apex  state  would  be  when 
displacing  the  actuator  in  the  following  time  step  and  in  the  second  half  of  the  stance,  the 
algorithm  is  essentially  creating  approximate  “grid  lines”  (see  Fig.  6.5  III)  that  are  based 
on  limited  numerical  simulations  of  the  stance  phase.  Computing  more  intermediate 
states,  i.e. ,  higher  order  fits  for  the  dynamics  of  the  system  subjected  to  actuation, 
would  certainly  result  in  a  better  estimate,  but  would  require  more  computational  time. 
For  applicability,  we  want  to  ensure  that  the  algorithm  finishes  well  within  the  time 
5t  between  actuation  updates.  As  a  result,  an  error  is  introduced  in  the  algorithm. 
Such  error  is  the  greatest  during  the  first  half  of  the  stance  phase,  when  we  have  more 
uncertainty  about  the  direction  the  motion  will  go,  and  results  in  the  bouncing  motion 
of  the  running  estimate  of  the  next  apex  state,  seen  in  Fig.  6.6  IV.  Fortunately,  the  first 
half  of  the  stance  is  also  characterized  by  a  higher  time-to-apex,  allowing  for  adequate 
time  to  significantly  reduce  the  approximation  error. 
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2.97  3.13  3..  28  3.44  3.6  2.94  3  3.07  3.13  3.2  3.26 

X  X 

Figure  6.5:  These  plots  show  an  example  of  the  steps  implemented  during  the  first 
half  of  the  stance  phase  in  order  to  approach  the  desired  apex  state. 
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Figure  6.6:  Plots  I  to  III  show  an  example  of  the  steps  implemented  during  the 
second  half  of  the  stance  phase  in  order  to  approach  the  desired  apex  state.  Plot  IV 
provides  a  summary  of  the  steps  implemented  during  the  entire  stance  phase. 
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6.2.1  Algorithm 

Let  us  call  P™ct(A tj)  and  —£™t(A tj)  the  maximum  actuator  movements  possible  in  the 
time  interval  A tj  with  respect  to  the  actuator’s  initial  position  and  its  specifications,  re¬ 
spectively  in  the  positive  (spring  compression)  and  negative  (spring  extension)  direction. 
Let  £®”td(A fj)  be  the  value  that  the  actuator  assumes  at  the  end  of  A tj. 

At  each  time  interval  A tn,  with  n  —  0,  1,  . . . ,  until  the  end  of  the  stance  phase: 

(i)  Compute  the  apex  state,  Sn+i,  that  the  SLIP  model  would  reach  if  keeping  the 
actuator  at  the  current  value  : 

Lctit)  =  t™?{A tn),  Vt  >  (n  +  1  )St, 

i.e. ,  the  apex  state  reachable  without  further  actuator  movements. 

•  If  £(t)  <  0  (first  half  of  stance  phase): 

(iiA)  Compute  the  apex  state  that  the  SLIP  model  would  reach  if  we  were  to  control 
the  actuator’s  movement  at  its  maximum  negative  or  positive  motion  for  a 
time  interval  At,  with  no  further  actuation  movement  occurring  for  the  rest 
of  the  stance  phase. 

±Ct(A*„+i),  Vt  e  Atn+1 

Lait)  =  < 

I  ^f(Atn+i),  otherwise. 

Let  us  call  these  points  S'rJ+1  and  S'“+1. 

(iiiA)  Compute  the  apex  state  that  the  SLIP  model  would  reach  if  we  were  to  move 

the  actuator  with  its  maximum  negative  or  positive  motion  for  a  time  interval 
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At  and  then  stop  until  the  beginning  of  the  second  half  of  the  stance  phase. 
Then,  the  actuator  moves  at  its  maximum  negative  or  positive  motion  until 
the  end  of  the  stance  phase. 


( 


4  d(t)  =  < 


±C*  (Atn+1), 
tn+1), 


eCct(A4+i), 


t  e  A tn+i 
Mt  :  l(t)  <  0 
Vf  :  £(t)  >  0, 


The  direction  of  movement  during  the  second  half  (£(f)  >  0)  is  dictated  by  the 
position  of  Sn+ 1  with  respect  of  Sdes,  as  shown  in  Fig.  6.7a.  In  fact,  during  the 
second  half  of  the  stance  phase,  a  positive  actuation  (+£4t)  would  increase 
both  y  and  x  values  at  the  next  apex,  while  a  negative  actuation  (— £™t)  will 
result  in  a  decrease  of  y  and  x.  Therefore,  if  Sn  is  in  zone  /  (see  Fig.  6.7a), 
then  £  =  1;  else,  £  =  —1.  Let  us  call  these  points  P„+i  and  Pn+V 

(ivA)  Consider  the  line  segments  s+  =  (5'rJ+1,  Pn+i)  and  s~  =  (S~+l,  Pn+l),  and 
compute  the  vector  distance  between  the  segments  and  the  desired  apex  state 
Sdes- 


d+  :=  vector  distance(s+,  Sdes ), 
d~  :=  vector  distance(s_,  Sdes )• 
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(vA)  Now,  Compute  the  next  actuator  value  desired,  £act(A £n+1),  as  follows. 


£act{Atn+l)  — 

|hCct(AW  1),  if  \d+\  <  \d~\ 

\-^Zt(Atn+ 1),  if|d“|<|d+| 

where  the  parameter  fi  is  a  weight  G  [0,  1].  In  particular,  if  the  two  vectors 
d+  and  dr  have  a  positive  dot  product,  then  n  —  1,  as  shown  in  Fig.  6.7c 
and  6.7d.  Otherwise,  the  desired  apex  point  is  ’’between’'  the  two  line  seg¬ 
ments,  and  the  required  actuator  movement  will  be  a  fraction  of  the  previous 
predictions,  as  shown  in  Fig.  6.7b: 


\d~ 

~  \d+ 1 

1  dr 

+  \d+ 1 
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Figure  6.7:  Subfig.  6.7a  shows  the  two  dividing  zones  where  Sn+i  can  be  with  respect 
to  the  position  of  Sdes ■  If  SVi+i  is  hr  zone  /,  during  the  second  half  of  stance  the  actu¬ 
ator  displacement  needs  to  increase;  else,  it  needs  to  decrease.  The  algorithm  creates 
approximate  grid  lines  (dotted  blue  line  segments)  based  on  limited  numerical  simula¬ 
tions  of  the  stance  phase:  only  the  states  S^+1)  S~+1,  P^+i:  and  Pn+i  are  computed. 
The  actuator  motion  during  the  first  half  of  stance  is  chosen  based  on  the  distance 
and  location  of  the  line  segments  with  respect  to  the  desired  state:  Subfig.  6.7b,  6.7c, 
and  6.7d  demonstrate  possible  location  of  vectors  d+  and  d~  with  respect  to  Sdes- 
Computing  more  states  would  result  in  a  better  fit  for  the  direction  along  which  the 
dynamics  of  the  system  move  when  subjected  to  actuation,  which  is  in  general  not 
a  straight  line,  and  as  a  consequence  it  would  result  in  a  more  accurate  choice  of 
actuator  motion.  However,  this  would  be  at  the  expense  of  the  computational  time, 
and  could  preclude  the  use  of  the  algorithm  in  real  time. 
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•  If  £(t)  >  0  (second  half  of  the  stance  phase): 

During  the  second  part  of  the  stance  phase,  both  x  and  y  functions  are  increasing 
or  decreasing,  as  seen  in  Fig.  6.2.  Therefore,  the  choice  of  whether  moving  the 
actuator  in  the  positive  or  negative  direction  is  again  done  based  on  the  position 
between  Sn+i  and  Sdes,  and  it  is  the  same  as  the  value  of  £  computed  in  (in)  . 


(iiB)  Compute  the  apex  state  that  the  SLIP  model  would  reach  if  we  were  to  control 
the  actuator’s  movement  at  its  maximum  velocity  for  the  time  interval  Atn+1 . 
Then,  no  further  actuation  movement  occurs  for  the  rest  of  the  stance  phase. 


Lct(t) 


UCct(^n+ 1),  VteAtn+1 
l£^f(Atn+1),  otherwise. 


Let  us  call  this  point  P^+i- 

(iiiB)  Compute  the  apex  state  that  the  SLIP  model  would  reach  if  we  were  to 
control  the  actuator’s  movement  at  half  of  its  maximum  (positive  or  negative) 
displacement  possible  for  a  time  interval  Atn+i  •  Then,  no  further  actuation 
movement  for  the  rest  of  the  stance  phase. 


Lct(t) 


Vt  G  Atn+i 
l^td(Atn+1),  otherwise. 


Let  us  call  this  point  Pn+V 

(ivB)  Consider  the  line  segments  s+  =  (Sn+ 1,  Pn+i)  and  s~  =  (Sn+ 1,  P~+  x ) .  Com¬ 
pute  the  vector  distance  between  the  segments  and  the  desired  apex  state 
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s. 


des • 


v+  :  =  vector  distance(s+,  Sdes), 
v~  :=  vector  distance(s_,  Sdes)- 


(vB)  Now,  compute  the  next  actuator  value  desired,  £act(Atn+i) ,  as  follows. 


^ act  (Afri.-|-l)  — 

\^act(Atn+ 1),  if|v+|<|u-| 

[/4Cct(AWi),  if|u“|<K| 

where  the  parameter  /i  is  a  weight  €  [0,  1].  In  particular,  if  both  or  neither 
vectors  are  orthogonal  to  the  respective  line  segment  s+  or  s~ ,  then  fi  —  1, 
as  shown  in  Fig.  6.8b  and  6.8c.  Else,  the  desired  apex  point  is  ’’between”  the 
two  line  segments,  and  the  required  actuator  movement  will  be  a  fraction  of 
the  previous  predictions,  as  in  Fig.  6.8a: 


V 

v+ 

v~ 

+ 

v+ 

(vi)  repeat  (i)  to  (v)  until  the  end  of  the  stance  phase. 

Note  that  the  choice  of  moving  the  actuator  with  its  maximum  velocity  and  accelera¬ 
tion  allowed  results  from  the  fact  that  a  higher  velocity  and  acceleration  correspond  to  a 
bigger  reachability  space  (e.g.,  see  Fig.  3.2  where  reachability  sets  computed  for  different 
values  of  maximum  velocities  are  compared).  The  algorithm  can  of  course  be  imple¬ 
mented  with  any  other  velocity  and  acceleration  values,  provided  they  are  in  accordance 
with  the  limitations  of  the  specific  actuator  used. 
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O  p+ 

4  rn+l 


P~  O  D+ 

n+10  W  rn+ 1 

4 


Figure  6.8:  These  figures  demonstrate  possible  locations  of  the  vectors  v+  and  v 
with  respect  of  Sdes- 
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6.3  Performance 

We  show  here  the  performances  of  our  algorithm.  In  particular,  we  start  with  con¬ 
siderations  on  the  required  time  delay  for  the  algorithm  implementation  and  the  time 
interval  used,  and  we  continue  with  simulation  results  showing  how  effective  the  algorithm 
is  in  reducing  the  distance  to  the  desire  state. 

6.3.1  Time  interval  and  delay  at  first  computational  instance 

The  choice  of  the  time  interval  5t  determines  the  frequency  at  which  the  algorithm 
is  updated.  Any  control  action  applied  in  real-time  during  the  stance  phase  needs  to 
be  executed  in  finite  time,  which  corresponds  to  the  duration  of  the  stance  phase  itself. 
The  execution  of  one  update  of  the  algorithm  is  only  a  function  of  the  processor  used, 
and  therefore  the  number  of  iterations  of  the  algorithm  depends  on  the  duration  of  the 
stance  phase,  which  in  turn  is  a  function  of  the  relative  spring  stiffness  7:  higher  values 
of  7  correspond  to  a  stiffer  spring,  and  therefore  a  shorter  stance  phase;  conversely,  lower 
values  of  7  correspond  to  a  longer  stance  phase  caused  by  a  softer  spring.  Let  us  choose 
initial  conditions  and  parameters  based  on  biological  data  for  a  typical  human,  as  in 
Table  6.1.  Then,  the  average  stance  time  is  shown  in  Fig.  6.9a  as  a  function  of  the 


Parameters 


9  = 

9.81  111/s2 

M  = 

80  kg 

4  = 

1  m 

(k,  0  = 

0.5  m 

(■act  ^ 

[-0.1,  0.1]  m 

y  e 

[1,  2.5]  m 

x  G 

[0.4,  5]  m/s 

Table  6.1:  Parameters  based  on  biological  data  for  a  typical  human, 
relative  spring  stiffness.  Considering  that  on  an  average  computer  (Intel  Core  i7  eight 
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core  processor  CPU,  2.8GHz)  running  Matlab  R2012a,  the  average  time  to  perform  all  the 
calculations  required  for  one  update  of  the  algorithm  is  ~  0.01,  0.02  [s],  one  can  ask  what 
values  of  7  for  a  real  robot  are  a  good  trade-off  between  number  of  algorithm  iterations 
and  considerations  on  biological  inspiration.  Fig.  6.9b  shows  how  many  algorithm  updates 
are  possible  on  average  during  stance,  for  5t  =  0.01  [s]  and  5t  =  0.02  [s] .  As  expected,  the 
smaller  the  value  of  7,  the  more  algorithm  iterations  can  be  performed.  A  study  by  [3] 
shows  that  the  relative  spring  stiffness  assumes  very  similar  values  for  various  gaits  in 
different  animals  (both  quadrupeds  and  bipeds).  For  runners,  values  were  found  between 
7  £  [7.1,  14.6],  while  7  £  [7.7,  13.6]  in  hoppers,  ft  is  then  reasonable  to  use  values  of 
7  £  [8,  13],  which  correspond  to  a  high  number  of  algorithm  iterations,  validating  the 
applicability  of  our  control  strategy. 
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Figure  6.9:  Subfig.  6.9a  shows  the  average  stance  time  as  a  function  of  7  for  a  set  of 
initial  apex  states  and  system’s  parameters  as  in  Table  6.1.  Subfig.  6.9b  shows  instead 
the  corresponding  number  of  algorithm  iterations. 
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When  implementing  a  control  action,  the  time  delay  of  such  controller  is  usually  a 
concern  that  requires  considerations.  In  the  case  of  the  algorithm  here  proposed,  each 
time  interval  A tn  is  dedicated  to  plan  the  desired  actuation  value  at  the  next  interval, 
£act(Atn+i) ■  This  implies  that  at  the  initial  time  interval  Ati  there  will  be  no  actuation; 
namely,  the  actuator  will  begin  its  motion  from  the  second  time-step.  We  want  here  to 
investigate  how  this  will  affect  the  performances  of  our  algorithm.  As  one  would  expect, 
increasing  time  delay  reduces  the  reachable  space  (Fig.  6.10)  by  limiting  the  maximum 
achievable  actuator  movement. 


Figure  6.10:  Reachability  sets  computed  for  different  delay  values.  As  expected, 
increasing  5t  reduces  the  reachable  space. 

The  delay  depends  on  the  time  interval  chosen  to  perform  our  algorithm,  St,  which  in 

turn  depends  on  the  computational  velocity  of  our  computer.  While  the  delay  may  seem 

a  limitation,  it  is  important  to  point  out  the  following.  On  an  average  computer,  the 

required  time  step  is  St  —  0.01  [s].  In  these  circumstances,  we  can  compare  the  reachable 

space  with  delay  with  the  reachable  spaces  obtained  with  other  control  strategies,  as 

108 


Reachability-based  control 


Chapter  6 


shown  in  Fig.  6.11.  As  we  can  see,  the  reachability  set  with  delay  is  still  bigger  (area- 
wise)  than  the  case  7 Zp  and  7ZS.  Thus  we  conclude  that  even  with  a  time  delay,  the 
proposed  algorithm  offers  significant  advantages  to  previous  control  strategies. 
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Figure  6.11:  Reachability  sets  with  5t  =  0.01  [s]  delay,  compared  to  7 Zp  and  7 Zs. 

This  performance  study  considers  the  ideal  SLIP  model,  where  the  body  is  a  point 
mass,  and  the  leg  is  massless.  Nonetheless,  we  expect  to  be  able  to  utilize  the  same 
control  strategy  on  a  real  platform  with  SLIP-like  dynamics.  Because  of  the  nature  of  its 
dynamics,  the  real  system  will  behave  similarly  to  the  ideal  model  when  subjected  to  a 
piston- like  actuation,  e.g.,  it  is  reasonable  to  believe  that  a  negative  actuation  will  result 
in  a  decrease  of  the  apex  height.  Additionally,  at  each  time  step  throughout  the  stance 
phase,  computations  of  the  subsequent  reachable  state  are  performed  numerically  using 
equations  of  motion  for  a  model  of  the  system.  Thus,  the  performance  of  the  algorithm 
is  strongly  tied  to  how  accurately  the  model  that  we  construct  of  the  physical  system 
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captures  its  actual  dynamics. 


6.3.2  Percentage  change 

We  will  evaluate  our  controller  using  equation  (6.4).  To  quantify  the  reduction  in 
Err ,  for  an  initial  apex  we  introduce  the  percentage  change,  PC ,  as: 


PC 


100- 


*S 'pass  Sdes  1 1 2  1 1  *-V  ^des  \  |  2 

1 1  Spass  *S ' des  \  |  2 


(6.5) 


where  Sdes  is  the  dimensionless  desired  apex  state,  Spass  is  the  dimensionless  state  reached 
without  actuation,  and  Sr  is  the  dimensionless  apex  state  reached  with  our  control  strat¬ 
egy.  Using  percentage  change  as  a  performance  index  allows  us  to  compare  the  system 
performance  with  and  without  control  action  over  a  broad  range  of  desired  states. 

The  performance  of  the  controller  has  been  tested  as  follows:  an  initial  apex  state  S'o 
and  touch-down  angle  9td  are  chosen,  and  a  set  of  10,  000  points  are  generated  inside  the 
corresponding  reachability  set  without  delay,  TZtot(So,  9td)-  Spass  is  defined  as  the  passive 
apex  state,  i.e.,  the  apex  state  that  is  reached  without  actuation  movement.  Then,  for 
each  point  chosen  inside  the  reachability  set,  the  controller  is  tested  by  setting  the  desired 
apex  state,  Sdes,  to  be  one  of  those  points,  with  a  controller  update  time  step  of  5t  =  0.01 
[s].  As  we  can  see  in  Fig.  6.12  and  Fig.  6.13,  the  effect  of  the  controller  is  to  dramatically 
reduce  the  error. 
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Figure  6.12:  Error  from  desired  point,  Err  =||  Sr—Sdes  || 2,  computed  for  10,  000  points 
randomly  taken  inside  the  reachable  set  TZ(So,  Otd},  with  So  =  {1.3  [m],  4  [m/s]}  and 
Otd  =  121.15  [deg].  The  update  time  has  been  chosen  to  be  5t  =  0.01  [s]. 
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Figure  6.13:  Percentage  change  graph  computed  for  10,000  points  randomly  taken 
inside  the  reachable  set  7 Z  =  {So,  0td}i  with  So  =  {1.3  [m],  4  [m/s]}  and  6td  =  121.15 
[deg].  The  update  time  is  St  =  0.01  [s],  and  the  inner  dotted  black  line  marks  the 
border  of  the  reachability  space  computed  with  a  St  gap.  The  black  x  mark  represents 
Spss,  the  apex  state  reached  without  any  actuator  motion.  As  expected,  the  percent 
change  of  the  subset  7 Z  —  7 Zgap  is  on  average  smaller  than  the  one  on  7 Z  alone,  since 
points  in  7Z  —  7Zgap  are  not  directly  reachable  with  our  proposed  strategy.  Note  that 
the  dark  (lower  improvement)  and  empty  areas  around  Spass  are  mainly  due  to  the 
fact  that  since  the  desired  apex  points  are  so  close  to  the  passively  reached  apex,  the 
margin  for  improvement  is  very  slim. 
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Figure  6.14:  Percentage  distribution  of  PC.  The  bars  represent  the  percentage  of 
points  with  a  PC  between  [0,  10],  (10,  20],...  (90,  100]  [%]  in  the  reachability  set 

without  delay  (light  grey  bars)  and  in  the  reachability  set  with  delay  St  =  0.01  [s] 
(dark  bars). 


Fig.  6.14  shows  how  many  point  inside  the  reachability  set  have  a  PC  value  between 
[0,  10],  (10,  20],...  (90,  100]  [%],  both  with  and  without  time  delay  (St  =  0.01  [s]), 

for  10,  000  random  desired  states.  If  we  consider  only  the  points  inside  the  time-delay 
reachability  set,  we  can  see  that  ~  50%  of  the  points  show  a  percentage  change  PC  G 
(90,  100],  meaning  that  the  actuator  motion  has  reduced  the  error  by  a  factor  between 
90%  and  100%.  However,  if  we  consider  instead  the  set  of  all  points  reachable  (i.e., 
without  delay),  the  percentage  of  points  with  PC  G  (90,  100]  is  instead  ~  42%,  as 
expected  smaller  than  its  subset  with  delay.  Indeed,  due  to  the  delay  in  the  execution  of 
the  controller,  points  outside  the  subset  71  —  7Zgap  are  not  reachable  anymore.  In  both 
cases,  though,  more  than  two  thirds  of  the  points  have  PC  G  (80,  100],  highlighting  the 
performance  of  our  proposed  controller,  even  if  it  requires  a  delay. 

Fig.  6.15,  and  Fig.  6.16a  and  6.16b  show  an  example  of  the  performance  of  our 
controller.  An  initial  apex  state  has  been  chosen  to  be  Sq  =  {1.3  [m],  4  [m/s]},  and 
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9td  =  121.15  [deg],  7  =  10.  The  desired  apex  state  is  Sdes  =  {1.5  [m],  3.8  [m/s]}. 
Maximum  velocity  and  acceleration  constant  were  chosen  to  be  vmax  =  0.5  [m/s]  and 
kacc  =  10  [m/s2].  After  running  the  controller  with  5t  =  0.01  [s],  the  reached  state  is 
Sr  =  {1.505  [m],  3.805  [m/s]}.  Fig.  6.15  shows  the  values  of  Sn  computed  at  each  time 
interval  At,  while  Fig.  6.16a  and  6.16b  show,  respectively,  the  errors  ||  Sdes  ~  Sn  H2  and 
its  x  and  y  component  at  each  time  interval,  and  the  actuator  value  £act(t)  over  time. 
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Figure  6.15:  Example  of  the  performance  of  the  proposed  control  algorithm.  The 
red  diamond  is  the  passive  state,  i.e.,  the  state  reached  without  any  actuation.  The 
x-marked  line  represents  the  apex  states  Sn,  n  =  1,2,...  computed  by  the  algorithm, 
while  the  yellow  circle  and  the  cyan  star  represent  respectively  the  desired  state,  Sdes, 
and  the  reached  state,  Sr. 
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Figure  6.16:  Subfig.  6.16a  and  6.16b  show  respectively  the  evolution  of  the  error  (6.4) 
and  the  evolution  of  its  component  in  x  and  y.  and  the  computed  actuation  for  the 
example  in  Fig.  6.15. 
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6.3.3  Percentage  change  over  a  set  of  consecutive  steps 

Section  6.3.2  presents  the  performance  of  our  control  strategy  for  one  specific  set  of 
initial  conditions:  starting  from  one  apex  state  and  the  touch-down  angle  necessary  to 
perform  a  symmetric  hop,  we  showed  the  performance  of  the  algorithm  when  trying  to 
reach  10,  000  states  inside  its  reachability  set.  To  show  that  the  presented  results  still 
hold  for  arbitrary  initial  conditions  and  touch-down  angles,  we  consider  a  path  of  1,  000 
subsequent  hops.  At  each  hop,  the  touch-down  angle  is  chosen  randomly  from  the  set  of 
angles  that  allow  a  successful  jump,  not  necessarily  symmetric,  where  a  successful  jump  is 
defined  as  any  apex-to-apex  hop  characterized  by  a  positive  velocity  and  an  apex  height 
exceeding  the  leg  length  to  avoid  collision  between  the  leg  and  the  ground.  Then,  the 
desired  apex  state  is  set  by  adding  to  each  dimension  of  the  passively  reached  state  an 
offset  picked  from  a  uniform  distribution.  Fig.  6.17  shows  all  the  apex  states  reached  by 
the  system  during  the  1,000  subsequent  hops.  As  each  one  of  these  states  is  the  initial 
condition  of  the  following  hop,  our  control  algorithm  has  been  tested  for  a  wide  variety 
of  initial  conditions  which  encompass  the  typical  operating  range  for  the  SLIP  model. 
Furthermore,  generating  desired  apex  states  by  adding  uniformly  distributed  noise  to  the 
passive  apex  state  ensures  that  our  controller  has  been  tested  for  a  wide  range  of  initial 
error.  Fig.  6.18  shows  the  resulting  percentage  change  for  1,000  hops.  As  we  can  see, 
over  two  thirds  of  the  states  have  a  percentage  change  PC  €  [90,100],  meaning  that 
the  controller  is  able  to  reduce  the  distance  to  the  desired  state  of  a  factor  between  90% 
and  100%.  Overall,  over  80%  of  the  states  have  a  PC  >  70.  More  in  detail,  Fig.  6.19 
shows  the  output  error  ||  Sr  —  Sdes  || 2  as  a  function  of  the  magnitude  of  the  initial  error 
I  Spass  —  Sdes  || 2-  For  each  initial  error  interval  of  size  0.025,  the  mean  value  (grey 
circle)  and  standard  deviation  (black  line)  of  the  output  error  has  been  computed:  as 
expected,  while  the  output  error  is  higher  for  higher  initial  error,  there  is  on  average  an 
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improvement  of  over  70%. 


Figure  6.17:  This  plot  shows  the  1,000  hops  (initial  conditions)  used  to  test  the 
controller.  The  initial  apex  state  has  been  chosen  as  So  =  {1-3  [m],  4  [m/s]}.  At  each 
hop,  the  desired  apex  state  has  been  chosen  by  adding  to  the  passively  reached  state  an 
offset  drawn  from  a  uniform  distribution  between  [—.3, ,  .3]  for  y  and  [^0.46,  0.46]  for 
x.  These  bounds  have  been  chosen  based  on  the  reachability  set  previously  computed 
for  the  initial  condition  So-  However,  note  that  each  hop  drives  the  system  to  a  new 
apex  state,  hence  a  new  reachability  set,  whose  computation  is  expensive  and  not 
feasible  in  real  time.  It  is  therefore  possible  at  any  given  hop  that  the  desired  apex  is 
not  in  the  corresponding  reachability  set. 
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Figure  6.18:  In  this  figure  the  resulting  percentage  change  values  are  shown:  as  we 
can  see,  675  points  have  a  PC  >  90,  i.e. ,  the  initial  distance  to  the  desired  apex 
is  reduced  of  more  than  90%.  Overall,  832  states  over  1,000  have  an  improvement 
above  70%.  Note  that  two  states  have  negative  percentage  change,  meaning  that  in 
these  two  cases  our  algorithm  failed  at  improving  the  initial  error,  due  to  numerical 
approximation  introduced  by  the  algorithm. 


'pass  & des  1 1 2 


Figure  6.19:  Plot  of  the  initial  error  versus  the  final  error.  Each  grey  circle  represents 
the  mean  value  of  the  final  error  for  all  the  initial  errors  in  the  corresponding  interval, 
i.e.,  between  [0, 0.025),  [0.025,  0.05),  etc.  The  black  vertical  lines  are  the  corresponding 
standard  deviations. 
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6.4  Conclusions 

In  this  chapter,  we  have  proposed  a  control  algorithm  for  the  actuated  SLIP  model 
to  reach  a  desired  apex  state.  Its  main  application  is  to  correct  errors  happening  at 
touch-down  and/or  during  the  stance  phase,  independent  of  any  foot  placement  strat¬ 
egy  employed.  Most  of  the  control  strategies  developed  in  the  past  for  the  SLIP  model 
(whether  passive  or  active)  are  aimed  at  controlling  the  touch-down  position  of  the  leg, 
9td ,  and/or  the  leg  length  at  touch-down,  in  order  to  reach  a  stable  gait  even  on  un¬ 
known  or  partially  known  terrain,  or  in  presence  of  system  noise.  For  the  actuated  case, 
this  choice  was  paired  with  actuation  displacement  during  the  stance  phase,  and  the 
desired  actuation  was  precomputed  during  flight  based  on  the  knowledge  on  the  terrain. 
Conversely,  our  strategy  focuses  on  reaching  a  desired  apex  state  without  any  insight 
on  the  choice  of  the  touch-down  angle:  the  actuation  is  computed  online  during  the 
stance  phase,  based  on  the  particular  state  of  the  system  at  each  time  interval.  The 
main  advantage  of  this  procedure  is  not  only  the  ability  to  counteract  sensing  errors  at 
landing  (e.g.,  ground  sensing  noise),  which  is  a  goal  taken  into  consideration  by  other 
works  as  well,  but  especially  the  ability  to  reduce  errors  and  disturbances  that  happen 
during  stance.  For  example,  a  slipping  of  the  foot,  or  an  external  force,  such  as  a  strong 
wind.  Furthermore,  we  want  to  point  out  that  our  strategy  does  not  aim  to  replace  any 
leg-placement  strategy:  on  the  contrary,  it  is  reasonable  to  think  that  our  controller  can 
easily  be  paired  with  any  other  leg-placement  or  path  planning  method. 

Another  contribution  of  this  chapter  is  a  qualitative  characterization  of  the  effects  of 
the  actuator  motion  on  the  reachable  apex  state.  Throughout  the  stance  phase,  comput¬ 
ing  the  spatial  relationship  between  the  current  reachable  state  without  further  actuation 
and  the  desired  state  gives  a  strong  indication  of  how  to  provide  additional  actuation  to 
the  system. 
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Finally,  by  delineating  the  reachability  space  for  a  given  apex  state  and  touch-down 
angle,  this  chapter  reveals  the  advantage  of  utilizing  the  actuation  at  its  maximum  capa¬ 
bilities.  In  fact,  we  show  that  updating  the  actuation  values  throughout  the  stance  phase 
allows  the  system  to  reach  a  wider  range  of  evenly  distributed  apex  states  with  respect 
to  the  ones  reached  by  pre-computing  some  fixed  actuator  motions. 
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Conclusions  and  future  work 


In  this  thesis  we  considered  the  SLIP  model,  used  as  a  template  for  modelling  the  dynam¬ 
ics  typical  of  legged  systems.  In  particular,  we  considered  an  actuated  version  of  the  SLIP, 
with  an  added  series  elastic  actuator  to  the  leg,  serving  the  purpose  of  adding/removing 
energy  to/from  the  system.  While  the  SLIP  model  has  been  a  topic  of  research  in  legged 
locomotion  for  several  decades,  studies  on  the  effect  of  actuation  on  the  system’s  behav¬ 
ior  are  still  not  complete.  Additionally,  most  of  the  proposed  control  strategies  focus  on 
steady-state  gaits  and  are  tuned  over  a  set  of  initial  conditions  or  terrain  profiles.  The 
purpose  of  this  thesis  is  to  fill  the  gap,  and  give  a  characterization  of  the  reachable  space 
for  the  active  SLIP,  paired  with  real-time  actuation  strategies  to  maximize  the  reachable 
space  and  drive  the  system  to  a  desired  state,  with  the  ability  to  change  stride  at  each 
step. 

Chapter  3  characterized  the  reachability  space  of  the  actuated  SLIP  model  by  acti¬ 
vating  the  series  elastic  actuator  at  any  possible  time  during  the  stance  phase.  Starting 
from  the  same  initial  conditions  at  apex,  for  each  possible  touch-down  angle  the  reach¬ 
able  space  is  shown  to  be  a  set  of  2-dimensional  manifolds,  one  per  touch-down  angle, 
forming  a  fan-shaped  3-dimensional  volume. 
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In  Chapter  4,  we  investigated  the  possibility  to  use  actuation  to  solve  in  closed-form 
the  system’s  dynamics.  This  has  been  done  by  enforcing  a  desired  trajectory  for  the 
dynamics  during  stance  phase.  The  actuator  only  acts  along  the  leg  of  the  system, 
and  therefore  it  is  not  possible  to  independently  enforce  a  trajectory  for  both  degrees 
of  freedom  (leg  angle  and  leg  length).  The  desired  trajectory  has  then  been  chosen  to 
be  consistent  with  the  system  dynamics.  In  particular,  the  desired  trajectory  has  been 
chosen  to  achieve  a  symmetric  motion,  i.e.,  apex  height  and  forward  velocity  are  preserved 
from  one  jump  to  the  next.  Enforcing  a  symmetric  trajectory  limits  the  reachable  space, 
which  becomes  a  point  in  the  3-dimensional  space.  Therefore  we  proposed  to  lock  the 
spring  during  stance  to  reach  an  asymmetric  take-off  state,  thus  combining  the  ability 
to  having  a  closed-form  solution  of  the  stance  phase  with  the  possibility  to  reach  several 
points  in  the  state  space. 

Chapter  5  aimed  to  extend  the  concept  of  utilizing  the  series  elastic  actuator  to  solve 
the  system’s  dynamics.  Via  partial  feedback  linearization  we  were  able  to  fully  solve  the 
dynamics  of  the  leg- length,  £(t),  and  we  then  proposed  an  approximation  of  the  leg-angle 
equation  that  took  into  account  the  actuator’s  dynamics.  The  main  benefit  of  having 
an  approximation  for  the  stance  phase  dynamics  is  the  reduced  computational  time  to 
predict  the  system’s  behavior  and  the  effect  that  diverse  actuation  strategies  have  on 
the  next  apex  reached.  We  then  proposed  a  two-element  control  strategy  to  apply  two 
actuation  references  throughout  the  stance  phase,  thus  driving  the  system  to  a  desired 
apex  state.  The  advantage  of  proposing  a  two-element  control  strategy  with  respect  to 
other  strategies  proposed  in  the  past  (for  example,  Raibert’s  thrust  at  mid  stance)  lays 
in  the  fact  that  it  is  possible  to  reach  a  wider  number  of  apex  states  in  one  step.  Also,  by 
having  an  approximation  of  the  stance  phase,  we  can  pre-compute  optimal  trajectories. 

Finally,  in  Chapter  6  we  studied  the  effect  that  applying  actuation  at  different  times 

during  stance  has  on  the  overall  reachable  region,  and  on  independently  controlling  each 
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of  the  two  apex  dimensions  (height  and  forward  velocity).  Based  on  an  understanding  of 
how  to  change  the  next  apex  through  actuation,  we  developed  an  algorithm  to  iteratively 
re-direct  the  dynamics  throughout  the  entire  stance  phase  to  achieve  a  desired  apex 
state,  if  reachable,  or  decrease  the  distance  to  it.  We  also  showed  how  this  strategy 
maximizes  the  reachable  set  with  respect  to  other  actuation  strategies,  and  since  it  is 
iteratively  computed  throughout  the  stance  phase,  it  can  deal  with  errors  or  adverse 
events  happening  during  motion. 

Overall,  the  main  contribution  of  this  thesis  is  to  understand  the  implications  that 
adding  a  series  elastic  actuator  have  on  the  system’s  dynamics,  and  particularly  on  the 
system’s  reachable  space,  and  to  determine  closed-form  solutions  for  the  stance  phase. 
Based  on  this  understanding,  we  proposed  real-time  control  strategies  for  actuator’s 
motion  that  drive  the  system  to  a  desired  state,  and  that  can  be  applied  to  a  variety  of 
situations,  from  running  on  rough  terrain  to  running  on  a  set  of  desired  footholds. 

Future  work  can  take  several  directions.  The  most  natural,  and  what  we  are  currently 
working  on,  is  the  application  of  our  proposed  strategies  to  a  hardware  implementation 
of  the  SLIP  model.  As  clearly  explained  in  [1],  the  ideal  SLIP  model  is  a  template 
that  can  then  be  ’’anchored”  to  systems  with  more  complex  dynamics.  The  prototype 
presently  under  development  and  study  at  the  UCSB  Robotics  Laboratory  can  be  seen 
in  Fig.  7.1a  and  7.1b.  The  UCSB  hopper  consists  of  a  body  of  mass  M  =  7.6  [kg] 
mounted  on  a  compliant  leg  of  natural  length  £0  =  0.54  [m]  and  mass  mieg  =  0.55  [kg] 
that  can  move  both  in  the  x  and  z  direction.  Actuation  on  the  leg  is  applied  through  a 
series  of  pulleys  activating  a  motor  at  the  hip  that  compresses/extends  the  spring.  As 
many  other  hopping/walking  robots,  the  UCSB  hopper  is  currently  mounted  on  a  boom 
to  reduce  lateral  stability  issues,  and  can  therefore  only  move  in  the  2-dimensional  space. 
However,  it  will  eventually  be  unmounted  from  the  boom  with  the  goal  to  control  it  in 
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(a) 


(b) 


Figure  7.1:  UCSB  hopper,  CAD  design  (a),  and  hardware  prototype  (b) 


the  3-dimensional  space.  The  body  is  free  to  rotate  around  the  hip,  where  a  motor  can 
apply  a  torque  to  control  leg  positioning  and  body  attitude.  However,  to  resemble  the 
point-mass  of  the  classical  SLIP  model,  the  body  is  currently  locked  via  a  mechanical 
stop  at  the  base  of  the  boom.  Encoders  at  the  leg,  actuators  and  boom  allow  us  to 
measure  leg’s  and  actuators’  positions.  A  gyroscope  on  board,  coupled  with  two  distance 
sensors  at  the  two  corners  of  the  body,  gives  us  a  measurement  of  the  body  tilt. 

There  are  several  challenges  in  applying  our  proposed  strategy  to  the  UCSB  hopper. 
First  of  all,  some  parameters  such  as  the  body  inertia  are  still  undetermined  and  need  to 
be  computed  via  system  identification.  Also,  the  mechanical  lock  is  applied  at  the  boom, 
thus  locking  the  body  on  the  boom  side.  However,  on  the  hopper  side  the  body  is  still  free 
to  vibrate  a  few  degrees,  both  due  to  impact  and  robot  movements  as  well  as  to  torsional 
dynamics  of  the  boom.  Sensor  noise  is  another  issue,  coupled  with  computing  the  exact 
instances  at  which  touch-down  and  take-off  happen.  The  control  strategies  proposed 
and  the  reachable  space  then  have  to  be  adapted  to  take  into  account  the  additional 
complexity  of  this  system  and  the  actuator’s  dynamics  versus  the  ideal  SLIP  model. 
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Other  directions  for  this  work  involve  extending  our  control  strategies  to  optimize 
energy  expenditure,  as  well  as  expanding  our  strategies  to  include  optimal  leg  placement 
at  touch-down.  It  would  be  interesting  as  future  work  to  quantify  the  range  of  terrain 
where  the  SLIP  model  successfully  works  (in  terms  of  metastability  and/or  energetic 
requirements).  For  example,  there  can  be  terrains  where  a  hopping  gait  is  optimal, 
but  on  other  terrain  profiles,  two  legs  and  a  walking  gait  would  be  more  efficient.  We 
can  then  foresee  the  possibility  and  the  benefits  of  combining  the  hopping  SLIP  model 
with  a  compass  gait  like  robot,  in  order  to  be  able  to  transition  between  walking  and 
running/hopping  and  negotiate  virtually  any  terrain  profile. 
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Appendix 

A.l  Planar  one-legged  hopper 

The  approximation  for  the  stance  phase  dynamics  discussed  in  Section  5.1  can  be 
extended  and  used  to  predict  the  dynamics  of  more  realistic  planar  one-legged  hoppers. 
In  particular,  we  show  here  how  our  approximation  can  be  applied  to  the  system  shown 
in  Fig.  5.19  and  described  in  Section  5.5.4. 

The  equations  of  motion  for  this  one-legged  hopper  can  be  easily  computed  from  its 
Lagrangian,  with  non  conservative  forces  given  by  the  torque  at  the  hip,  ThiP-  During 
flight,  the  center  of  mass  of  the  system  follows  a  ballistic  trajectory,  and  the  servo  at  the 
hip  moves  the  body  to  drive  the  leg  to  a  desired  touch-down  angle.  The  initial  conditions 
for  the  leg-angle  and  leg-length  dynamics  during  the  stance  phase  depend  on  the  state 
of  the  body  at  touch-down,  Xbodyitro),  VbodyitTD),  and  their  velocities.  These  values  can 
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be  approximated  from  the  touch-down  state  of  the  center  of  mass  as: 


Xbody(tTD)  ~  %com  {tm)  + 
UbodyitTD )  ~  Vcom^TD)  + 


m(£o  -  £i) 

M  +  m 
m(£0  -  di) 
M  +  m 


sin  0TD, 

COS 6TDi 


Xbody(tTD )  ~  i'com  (tm), 


j)body{tTD )  ~  Vcom 


The  energy  loss  at  impact  is  taken  into  account  at  take-off  as  follows: 

.  .  M.  m£i  .  na 

XcomV-To)  ~  - xbody(tTO)  +  — — - sm(6hx>)nro; 

M  +  m  M  +  m 

.  M.  mil  na 

VcomitTO)  ~  T-J— ; - ybodyV-TO)  +  "7T  J - COS  ( 0TO)0TO ■ 

M  +  m  M  +  m 

During  stance,  the  partial  feedback  linearization  (5.1)  is  applied  to  the  system;  hence, 
the  closed  form  solution  for  the  leg  length  i  derived  in  section  5.1  is  still  valid.  The  series 
elastic  actuator  as  a  function  of  input  current  ucurr(t )  is  modelled  as: 

La  =  ~  (  b£act  -  Ku(t))  -  fsign(Lct),  (A.l) 

meff 

where  me//  is  the  effective  mass  seen  by  the  actuator,  b  is  the  damping  constant,  /  is 
the  friction  constant,  and  ft  is  the  motor  constant. 

In  general,  the  non-zero  momentum  of  the  body  affects  the  dynamics  of  the  leg-angle, 
6(t).  Then,  9  can  be  approximated  with  the  solution  computed  in  section  5.1,  6appr ,  with 
the  addition  of  the  term  that  includes  the  effect  of  the  body  motion  over  the  leg  angle, 

Of 

9{t)  =  9appr(t)  +  (t),  (A. 2) 
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where 

n  J(M  +  m)  v 

0  Ji(M  +  m)  +  Mm (£  -  4  )■ 2 

From  the  equations  of  motion,  we  have  that  Jcft  =  Tj-np\  hence,  by  choosing  an  appropriate 
torque  at  the  hip,  6 $  can  be  analytically  solved.  In  particular,  to  avoid  excessive  tipping 
of  the  body,  tup  has  been  chosen  as  the  constant  value  necessary  to  keep  the  body 
upright.  The  solution  for  9(f)  computed  in  (A. 2)  is  an  approximation.  During  stance, 
a  high  level  PID  controller  for  ThiP  is  applied  to  control  the  system  to  follow  the  target 
dynamics  (A. 2). 


Simulation  Parameters 


M  = 

10  kg 

m  = 

1  kg 

J  = 

10  kg  •  m2 

Ji  = 

1  kg  •  m2 

4  = 

1  m 

h  = 

0.5  m 

Pact  C 

[-0.1,  0.1]  m 

4  e 

[—0.05,  0.05]  m 

Vact  — 

1  m/s 

Vc  = 

0.5  m/s 

7  = 

20 

^ curr  ^ 

[-20,  20]  A 

7~hip  ^ 

[—50,  50]  N  •  m 

b  = 

50 

f  = 

5 

K  = 

68.8 

meff  = 

10  kg 

Table  A.l: 


This  approximation  has  been  tested  on  a  path  of  100  jumps  with  varying  target  apex 

states,  as  described  in  5.5.4,  with  parameters  from  Table  A.l.  Fig.  A. la  show  an  example 

of  the  approximated  dynamics  of  the  leg  angle  during  stance  compared  to  their  numerical 

solution,  and  A. lb  the  hip  torque  required.  Fig.  A. 2  shows  the  series  elastic  actuator 

motion,  £act  and  its  two  components:  ini  and  £Cl  and  the  current  required  ucurr. 
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(a) 


(b) 

Figure  A.l:  (a)  Data  showing  the  approximated  solution  of  the  dynamics  of  Q(t)  and 
9(t)  compared  to  the  numerically  computed  solution.  Fig.  (b)  shows  the  required 
torque  at  the  hip. 


In  Fig.  A. 3a  we  can  see  the  approximation  of  the  trajectory  of  the  body  mass  com¬ 
pared  to  the  numerical  computed  solution,  while  Fig.  A. 3b  show  the  approximation  and 
numerical  solution  of  the  forward  velocity  of  the  center  of  mass,  for  a  few  jumps. 
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Figure  A. 2:  Series  elastic  actuator  motion  £act(t )  =  £ni  +  lc  and  current  ucurr.  A 
PID  controller  drives  the  actuator  dynamics  (A.l)  to  the  target  values.  The  constant 
values  desired  during  the  first  and  second  half  of  the  stance  phase  are,  respectively, 
£c\  =  0.03,  and  =  0.01. 
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(a) 


(b) 

Figure  A. 3:  Approximation  of  the  trajectory  (a)  and  forward  velocity  ( b )  of  the  mass 
compared  to  the  numerically  computed  solution.  The  apparent  mismatch  during  the 
stance  phase  in  x  is  due  to  the  fact  that  the  energy  loss  at  impact  is  modelled  at 
take-off  in  the  approximated  solution. 


131 


Bibliography 


[1]  R.  J.  Full  and  D.  E.  Koditschek,  Templates  and  anchors:  Neuromechanical 
hypotheses  of  legged  locomotion ,  Journal  of  Experimental  Biology  202  (1999) 
3325-3332. 

[2]  R.  Blickhan,  The  spring-mass  model  for  running  and  hopping ,  Journal  of 
Biomechanics  22  (1989),  no.  11/12  1217-1227. 

[3]  R.  Blickhan  and  R.  J.  Full,  Similarity  in  multilegged  locomotion:  Bouncing  like  a 
monopode ,  Journal  of  Comparative  Physiology  A:  Neurothology,  Sensory,  Neural, 
and  Behavioral  Physiology  173  (Nov.,  1993)  509-517. 

[4]  J.  Schmitt  and  P.  Holmes,  Mechanical  models  for  insect  locomotion:  dynamics  and 
stability  in  the  horizontal  plane  i.  theory ,  Biological  Cybernetics  83  (2000),  no.  6 
501-515. 

[5]  J.  Schmitt  and  P.  Holmes,  Mechanical  models  for  insect  locomotion:  dynamics  and 
stability  in  the  horizontal  plane  ii.  application ,  Biological  Cybernetics  83  (2000), 
no.  6  517-527. 

[6]  R.  Full,  C.  Farley,  and  J.  Winters,  Musculoskeletal  dynamics  in  rhythmic  systems: 
A  comparative  approach  to  legged  locomotion ,  in  Biomechanics  and  Neural  Control 
of  Posture  and  Movement  (J.  Winters  and  P.  Crago,  eds.),  pp.  192-205.  Springer 
New  York,  2000. 

[7]  M.  A.  Daley,  G.  Felix,  and  A.  A.  Biewener,  Running  stability  is  enhanced  by  a 
proximo- distal  gradient  in  jointy  neuromechanical  control ,  Journal  of  Experimental 
Biology  210  (2007)  383-394. 

[8]  D.  Dudek  and  R.  J.  Full,  Passive  mechanical  properties  of  legs  from  running 
insects ,  Journal  of  Experimental  Biology  209  (2006)  1502-1515. 

[9]  S.  Sponberg  and  R.  J.  Full,  Neuromechanical  response  of  musculo -skeletal 
structures  in  cockroaches  during  rapid  running  on  rough  terrain ,  Journal  of 
Experimental  Biology  211  (2008)  433-446. 

[10]  M.  H.  Raibert,  Legged  Robots  that  Balance.  MIT  Press,  1986. 


132 


[11]  G.  Zeglin,  The  Bow  Leg  Hopping  Robot.  PhD  thesis,  Carnegie  Mellon  University, 
Pittburgh,  PA,  USA,  Oct.,  1999. 

[12]  M.  Ahmadi  and  M.  Buehlcr,  Stable  control  of  a  simulated  one-legged  running  robot 
with  a  hip  and  leg  compliance ,  IEEE  Transactions  on  Robotics  13  (1997),  no.  1. 

[13]  J.  Schmitt  and  J.  Clark,  Modeling  posture- dependent  leg  actuation  in  sagittal  plane 
locomotion ,  Bioinspiration  and  Biomimetics  4  (2009)  1-17. 

[14]  B.  Andrews,  B.  Miller,  J.  Schmitt,  and  J.  E.  Clark,  Running  over  unknown  rough 
terrain  with  a  one-legged  planar  robot,  Bioinspiration  &  Biomimetics  6  (2011), 
no.  2  026009. 

[15]  M.  Rutschmann,  B.  Satzinger,  M.  Byl,  and  K.  Byl,  Nonlinear  model  predictive 
control  for  rough  terrain  hopping ,  in  IEEE/RSJ  Int.  Conf.  on  Intelligent  Robots  & 
Systems,  (Villamoura,  Algarve,  Portugal),  pp.  1859-1864,  Oct.,  2012. 

[16]  K.  Byl,  M.  Byl,  M.  Rutschmann,  B.  Satzinger,  L.  van  Blarigan,  G.  Piovan,  and 
J.  Cortell,  Series- elastic  actuation  prototype  for  rough  terrain  hopping,  in  IEEE 
International  Conference  on  Technologies  for  Practical  Robot  Applications, 

pp.  103-110,  2012. 

[17]  J.  Seipel  and  P.  Holmes,  A  simple  model  for  clock- actuated  legged  locomotion, 
Regular  and  Chaotic  Dynamics  12  (2007),  no.  5  502-520. 

[18]  H.  R.  Vejdani  and  J.  W.  Hurst,  Swing  leg  control  for  actuated  spring-mass  robots, 
in  Int.  Conf.  on  Climbing  and  Walking  Robots  and  the  Support  Technologies  for 
Mobile  Machines ,  (Baltimore,  MD,  USA),  pp.  536-542,  2012. 

[19]  S.  Riese  and  A.  Seyfarth,  Stance  leg  control:  variation  of  leg  parameters  supports 
stable  hopping,  Bioinspiration  &  Biomimetics  7  (2012),  no.  1  016006. 

[20]  A.  Seyfarth,  H.  Geyer,  M.  Gunther,  and  R.  Blickhan,  A  movement  criterion  for 
running,  Journal  of  Biomechanics  35  (2002)  649-655. 

[21]  R.  M.  Ghigliazza,  R.  Altendorfer,  P.  Holmes,  and  D.  Koditschek,  A  simply 
stabilized  running  model,  SIAM  Rev.  47  (Mar.,  2005). 

[22]  A.  Seyfarth,  H.  Geyer,  and  H.  Herr,  Swing-leg  retraction:  a  simple  control  model 
for  stable  running,  Journal  of  Experimental  Biology  206  (2003)  2547-2555. 

[23]  J.  G.  D.  Karssen,  M.  Haberland,  M.  Wisse,  and  S.  Kim,  The  optimal  swing-leg 
retraction  rate  for  running,  in  IEEE  Int.  Conf.  on  Robotics  and  Automation, 
(Shangai,  China),  pp.  4000-4006,  May,  2011. 


133 


[24]  R.  Altendorfer,  D.  E.  Koditschek,  and  P.  Holmes,  Stability  analysis  of  legged 
locomotion  models  by  symmetry-factored  return  maps,  International  Journal  of 
Robotics  Research  23  (2004),  no.  10-11  979-999. 

[25]  M.  Ernst,  H.  Geyer,  and  R.  Blickhan,  Spring-legged  locomotion  on  uneven  ground: 
a  control  approach  to  keep  the  running  speed  constant,  in  12th  International 
Conference  on  Climbing  and  Walking  Robots  (CLAWAR),  pp.  639-644,  2009. 

[26]  M.  Ahmadi  and  M.  Buehler,  Controlled  passive  dynamic  running  experiments  with 
the  ARL-Monopod  II,  IEEE  Transactions  on  Robotics  22  (2006),  no.  5. 

[27]  J.  K.  Hodgins  and  M.  H.  Raibert,  Adjusting  step  length  for  rough  terrain 
locomotion,  IEEE  Transactions  on  Robotics  and  Automation  7  (1991),  no.  3. 

[28]  D.  Koepl  and  J.  Hurst,  Impulse  control  for  planar  spring-mass  running,  Journal  of 
Intelligent  and  Robotic  Systems  74  (2014)  589-603. 

[29]  O.  Arslan  and  U.  Saranli,  Reactive  planning  and  control  of  planar  spring-mass 
running  on  rough  terrain,  IEEE  Transactions  on  Robotics  28  (2012),  no.  3  567-579. 

[30]  E.  T.  Whittacker,  A  treatise  on  the  analytical  dynamics  of  particles  and  rigid 
bodies.  Cambridge  University  Press,  New  York,  fourth  ed.,  1904. 

[31]  W.  J.  Schwind  and  D.  E.  Koditschek,  Approximating  the  stance  map  of  a  2  DOF 
monoped  runner,  Journal  of  Nonlinear  Science  10  (2000),  no.  5  533-588. 

[32]  H.  Geyer,  A.  Seyfarth,  and  R.  Blickhan,  Spring-mass  running:  simple  approximate 
solution  and  application  to  gait  stability,  Journal  of  Theoretical  Biology  232  (Feb., 
2005)  315-328. 

[33]  O.  Arslan,  U.  Saranli,  and  O.  Morgiil,  Approximate  stance  map  of  the  spring  mass 
hopper  with  gravity  Correction  for  Nonsymmetric  Locomotions,  in  IEEE  Int.  Conf. 
on  Robotics  and  Automation,  (Kobe,  Japan),  pp.  2388-2393,  May,  2009. 

[34]  H.  Yu,  M.  Li,  and  H.  Cai,  Approximating  the  stance  map  of  the  SLIP  runner  based 
on  perturbation  approach,  in  IEEE  Int.  Conf.  on  Robotics  and  Automation,  (Saint 
Paul,  MN,  USA),  pp.  4197-4203,  May,  2012. 

[35]  Z.  Shen  and  J.  Seipel,  A  piecewise-linear  approximation  of  the  canonical 
spring-loaded-inverted-pendulum  model  of  legged  locomotion,  Journal  of 
Computational  and  Nonlinear  Dynamics  (January  2015). 

[36]  J.-C.  Zufferey,  First  jumps  of  the  3d  bow  leg  hopper,  Master’s  thesis,  E’cole 
Polytechnique  Fe’de’rale  de  Lausanne,  and  Carnegie  Mellon  University,  2001. 


134 


[37]  G.  Zeglin  and  H.  B.  Jr.,  First  hops  of  the  3d  bow  leg,  in  Int.  Conf.  on  Climbing  and 
Walking  Robots  and  the  Support  Technologies  for  Mobile  Machines,  pp.  357-364, 
2002. 

[38]  U.  Saranli,  M.  Buehler,  and  D.  E.  Koditschek,  RHex:  A  simple  and  highly  mobile 
hexapod  robot,  International  Journal  of  Robotics  Research  20  (2001)  616-631. 

[39]  P.  Gregorio,  M.  Ahmadi,  and  M.  Buehler,  Design,  control,  and  energetics  of  an 
electrically  actuated  legged  robot,  Systems,  Man,  and  Cybernetics,  Part  B: 
Cybernetics,  IEEE  Transactions  on  27  (Aug,  1997)  626-634. 

[40]  J.  W.  Hurst,  J.  Chestnutt,  and  A.  Rizzi,  Design  and  philosophy  of  the  BiMASC,  a 
highly  dynamic  biped,  in  IEEE  Int.  Conf.  on  Robotics  and  Automation,  (Roma, 
Italy),  Apr.,  2007. 

[41]  B.  Dadashzadeh,  H.  R.  Vejdani,  and  J.  Hurst,  From  template  to  anchor:  A  novel 
control  strategy  for  spring-mass  running  of  bipedal  robots,  in  IEEE/RSJ  Int.  Conf. 
on  Intelligent  Robots  &  Systems,  (Chicago,  IL,  USA),  pp.  2566-2571,  Sept.,  2014. 

[42]  T.  Lejeune,  P.  Willems,  and  N.  Heglund,  Mechanics  and  energetics  of  human 
locomotion  on  sand,  Journal  of  Experimental  Biology  201  (1998)  2071-2078. 

[43]  O.  Arslan,  U.  Saranli,  and  O.  Morgiil,  Reactive  footstep  planning  for  a  planar 
sprinq  mass  hopper,  in  IEEE/RSJ  Int.  Conf.  on  Intelligent  Robots  &  Systems ,  (St. 
Louis,  USA),  Oct.,  2009. 

[44]  G.  Piovan  and  K.  Byl,  Enforced  symmetry  of  the  stance  phase  for  the  spring-loaded 
inverted  pendulum,  in  IEEE  Int.  Conf.  on  Robotics  and  Automation,  (Saint  Paul, 
MN,  USA),  pp.  1908-1914,  May,  2012. 


135 


